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FOREWORD 


The  Model  Output  Statistics  (MOS)  approach  has  been  used  extensively  to  pre¬ 
dict  a  variety  of  meteorological  phenomena.  The  approach  consists  of  extract¬ 
ing  statistical  relationships  between  a  number  of  model-generated  variables 
(the  predictors)  and  the  forecast  variable  of  interest  (the  response  —  in  our  case, 
cloud  cover  fraction).  The  model  data  are  typically  output  frum  a  physically  based 
(i.e.,  numerical)  forecast  model  such  as  those  used  at  the  U.S.  National  Weath¬ 
er  l^rvice. 

In  order  to  develop  a  MOS-based  forecasting  procedure  with  sufficient  skill  to 
be  of  value  to  users,  a  tremendous  amount  of  data  must  be  acquired  and  ana¬ 
lyzed.  Furthermore,  those  data  must  be  collected  from  the  same  physical  envi¬ 
ronment  in  which  the  model  vdll  be  applied.  Thus,  a  large  collection  of  winter 
data  will,  almost  surely,  be  of  little  value  for  forecasting  weather  variables  dur¬ 
ing  the  summer.  Similarly,  if  the  large  physical  model  whose  outputs  include 
some  of  the  predictors  is  changed  for  any  reason,  then  MOS  relationships  must 
be  re-derived,  or  at  least  adjusted  and  verified,  to  accommodate  the  changes. 

In  spite  of  the  above  complexities,  the  possibility  of  large-scale  applications  of 
the  MOS  approach  to  cloud-cover  fraction  forecasting  is  attractive  to  the  U.S. 
Air  Force  broause  of  its  extreme  efficiency  once  the  MOS  functional  relation¬ 
ship  has  been  established.  This  project,  therefore,  was  funded  to  perform  a 
carefully  designed  feasibility  study  of  the  MOS  approach  to  cloud-cover  frac¬ 
tion  forecasting.  The  study  considered  forecasts  of  two  response  variables: 
cloud-cover  fraction  and  probability  of  cloud-cover  fraction  less  than  a  specified 
threshold  value.  The  latter  response  variable  is  of  special  interest  to  one  part 
of  the  customer  community. 

The  one  year  of  data  made  available  to  TASC  by  the  government  for  this  feasibil¬ 
ity  study  would  not  be  adequate  to  support  frill  development  of  an  operational 
MOS  forecasting  system  valid  over  future  years.  It  was,  however,  adequate  to 
judge  feasibility  in  terms  of  potential  accuracy  of  the  forecasts,  variability  of  the 
forecast  equations  with  seasons  and  with  physical  locations,  and  similar  issues. 
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INTRODUCTION 


1.1  OVERVIEW  OF  OUR  GOALS  AND  METHODOLOGY 

The  goal  of  this  project  is  to  “demonstrate  whether  or  not  the  model  output  statistics 
(MOS)  approach  should  be  considered  as  a  serious  contender  for  a  low-cost,  low-risk  im¬ 
provement  to  the  current  Air  Force  Global  Yieather  Central  (GWC)  cloud  amount  forecasting 
system”  (Ref.  26).  The  Department  of  Defense  has  assigned  GWC  the  task  of  providing  accu¬ 
rate,  timely,  global  cloud  forecasts  at  three-hour  intervals  out  to  48  hours  (Ref.  9).  Reduced 
visibility  due  to  clouds  in  the  fields  of  view  of  electro-optical  sensor  systems  adversely  im¬ 
pacts  military  operations  such  as  reconnaissance  missions,  tactical  air  strikes,  and  mid-air 
refueling.  lb  maximize  the  chance  of  mission  success.  Army  and  Air  Force  decision  makers 
require  accurate,  timely  forecasts  of  cloud  amount.  In  addition,  mission  planning  may  re¬ 
quire  forecasts  of  the  probability  of  exceeding  some  mission-specific  cloud  coverage. 

During  this  project,  TASC  developed  model  output  statistics  regression  models  for 
3-,  6-,  and  9-hour  cloud  amount  forecasts  and  compared  the  results  to  persistence  fore¬ 
casts  to  determine  the  utility  of  the  MOS  approach  for  short-range  cloud  prediction.  Re¬ 
gression  models  were  developed  for  up  to  eight  different  times  of  day  on  the  eighth-mesh 
grid  (approximately  40  km  grid  spacing  at  midlatitudes)  within  climatol<^cally  distinct 
regions  of  the  northern  hemisphere  during  all  seasons.  As  data  sources  for  the  MOS  fore¬ 
casts,  we  used  total  cloud  amount  fields  fi^m  the  Real-Time  NEPHanalysis  (RTNEPH) 
model  and  meteorological  variables  from  the  Air  Force  Global  Spectral  Model  (GSM).  We 
also  used  terrain  type  and  elevation  variables  from  a  terrain  database  (on  the  eighth- 
mesh  grid)  compiled  at  Phillips  Laboratory.  The  RTNEPH  cloud  amount  fields  were  also 
used  as  "ground  truth.” 

The  modeling  process  consisted  of  the  following  major  components.  A  more  com¬ 
plete  description  of  the  process  can  be  found  in  Section  3. 

STEP  1  Construct  clusters  of  RTNEPH  grid  points  based  on  cloud  climatology 
and  data  availability.  Select  clusters  for  which  to  develop  regression  models. 

STEP  2  Develop  a  pool  of  potential  predictors  from  RTNEPH  and  GSM  model 
data.  This  pool  includes  raw  and  derived  variables  along  with  averaged,  differ¬ 
enced,  and  transformed  variables. 

STEP  3  Develop  regression  models  for  total  cloud  amount  using  variables  se¬ 
lected  from  the  pool  of  predictoro.  Determine  conditional  probabilities  of  cloud 
amount  thresholds. 
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STEP  4  Predict  MOS  model  performance  for  future  data  sets  using  the  jack¬ 
knife  resampling  technique  to  analyze  total  cloud  cover  forecasts  and  a  Bayesian 
approach  to  analyze  probabilistic  forecasts. 

Although  the  overall  process  is  shown  as  a  linear  sequence,  we  continued  to  iterate 
on  the  sequence  as  the  modeling  effort  proceeded.  Keeping  in  mind  the  goal  of  the  project 
— a  feasibility  study — we  emphasized  research  and  exploration  and  a  d3mamic  modeling 
process  in  which  findings  from  later  steps  led  us  to  modify  earlier  ones. 

Although  not  tasked  to  develop  exhaustive  operational  MOS-based  forecast  models, 
wb  were  acutely  aware  of  the  operational  constraints  of  the  RTNEPH  and  GSM  model  runs. 
Therefore  we  were  careful  to  include  only  those  data  sources  that  are  operationally  available 
at  any  given  time  as  MOS  inputs.  Our  development  data  set  consisted  of  one  year  (1989)  of 
Global  Spectral  Model  output  mapped  to  RTNEPH  coordinates  within  20  RTNEPH  boxes  in 
the  northern  hemisphere  and  the  RTNEPH  doud  analysis  grids  for  the  same  time  period  and 
spatial  regions.  Both  of  these  data  sets  are  described  further  in  Section  2.  Jacks  et  al.,  in 
Ref.l7,  point  out  that  at  least  two  years  of  model  data  are  needed  to  develop  meaningful  MOS 
equations  based  on  the  Nested  Grid  Model.  We  use  one  year  of  data  in  this  feasibility  study 
with  the  understanding  that  more  robust  models  that  better  account  for  year-to-year  meteo¬ 
rological  variability  would  be  necessary  in  an  operational  setting. 


1.2  OVERVIEW  OF  THIS  REPORT 

The  remainder  of  this  section  discusses  the  current  forecasting  capabilities  at  GWC 
and  the  potential  improvements  to  be  made  by  using  a  MOS  approach.  Section  2  provides 
descriptions  of  the  RTNEPH  and  GSM  data  sources  used  in  model  development.  Section  3 
contains  a  description  of  our  modeling  process  along  with  the  modeling  and  validation  re¬ 
sults.  Finally,  we  summarize  and  make  recommendations  for  further  research  in  Section  4. 


1.3  THE  CURRENT  GWC  CLOUD  FORECASTING  SYSTEM 

At  present,  GWC  operates  a  Global  Spectral  numerical  prediction  Model  (GSM),  a 
global  cloud  analysis  model  (RTNEPH),  a  global  High  Resolution  Analysis  System  (HIRAS) 
to  initialize  the  GSM,  and  three  cloud  forecasting  models:  the  IVopical  Cloud  Forecast  Mod¬ 
el  (TRONEW),  the  Five-Layer  Cloud  Forecast  Model  (SLAYER),  and  the  H^h-Resolution 
Cloud  Prognosis  Model  (HRCP).  Figure  1  depicts  the  system’s  component  models  and  the 
data  flow  among  them.  The  RTNEPH  cloud  analysis  model  provides  initial  data  for  the 
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Figure  1  Global  Weather  Central  Cloud  Analysis  and  Forecast  System 

three  cloud  forecast  models.  It  and  the  Global  Spectral  Model  are  discussed  in  some  detail 
in  Section  2.  Here  we  summarize  the  three  cloud  forecast  models. 

IVopical  Cloud  Forecast  Model  (TRONEW) — ^The  TRONEW  model  (Ref.  9)  is 
based  on  RTNEPH  cloud  analyses.  It  persists  the  ditimal  cycle  for  24  hours,  and  outputs 
total  cloud  amounts  and  layer  cloud  amounts  and  types  at  three-hour  intervals  out  to  21 
hours.  The  complete  half-mesh  grid  (approximately  100  km  grid  spacing  at  the  equator) 
is  used  in  both  hemispheres  to  provide  global  coverage.  Known  problems  include:  the  sim¬ 
ple  forecast  technique  (no  dynamics),  poor  horizontal  resolution,  and  poor  temporal  reso¬ 
lution  in  the  tropics  due  to  inadequate  DMSP  satellite  coverage  in  that  area. 

Five-Layer  Cloud  Forecast  Model  (SLAYER) — ^The  SLAYER  model  (Ref.  9)  in¬ 
puts  are  the  RTNEPH  cloud  analyses  and  the  HIRAS/GSM  analyses  and  forecasts  of 
wind,  temperature,  and  moisture.  It  uses  a  quasi-Lagrangian  trajectory  scheme  to  advect 
temperature  and  condensation  pressure  spread  (CPS).  A  subset  of  the  half-mesh  grid  (an 
octagon  which  excludes  latitudes  between  0  and  approximately  13  degrees)  is  used  in  each 
hemisphere.  Forecasts  are  made  on  five  vertical  levels:  gradient  (60  mb  above  the  sur¬ 
face),  850, 700, 500,  and  300  mb.  Output  is  generated  every  three  hours  to  48  hours  in  the 
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northern  hemisphere  and  to  24  hours  in  the  southern  hemisphere.  Output  parameters  in¬ 
clude:  total  cloud  amounts,  layered  cloud  amounts,  precipitation  amounts  and  types,  ac¬ 
cumulated  precipitation,  icing  type,  dew  point  depression,  trajectories,  Showalter 
stability,  temperature,  and  CPS.  Known  problems  include:  lack  of  coverage  near  the  equa¬ 
tor,  poor  performance  in  the  boundaiy  layer,  lack  of  entrainment  and  cumulus  develop¬ 
ment,  poor  handling  of  cirrus  and  fog/stratus  discrimination,  and  the  simplicity  of  the 
empirical  relationship  between  CPS  and  cloud  amount. 

High-Resolution  Cloud  Prognosis  Model  (HRCP)— The  HRCP  model  (Ref.  9) 
inputs  are  the  RTNEPH  analyses,  the  SLAYER  CPS  and  trajectory  fields,  and  terrain 
heights.  It  can  use  either  persistence  or  forecast  trends  in  CPS  to  project  CPS  out  to  three, 
six,  and  nine  hours.  CPS  is  then  related  to  cloud  amount  using  empirical  curves  originally 
developed  in  the  1960’s.  The  HRCP  model  operates  in  a  window  on  the  eighth-mesh  grid. 
The  window  is  movable  in  the  sense  that  up  to  a  maximum  of  13  of  the  120  RTNEPH  boxes 
(60  per  hemisphere)  may  be  specified  for  a  given  run.  Each  box  contains  4096  (64  by  64) 
points  with  a  nominal  spacing  of  40  km.  Outputs  include  total  cloud  amounts  and  cloud 
amounts  and  types  for  low,  middle,  and  high  layers.  Known  problems  include:  low  quality 
forecasts  for  the  boundary  layer  and  the  middle  layers,  and  inadequate  parameterization 
of  cloud  amount  as  a  function  of  CPS. 

1.4  STATISTICAL  CLOUD  FORECASTING 

1.4.1  Terminology 

The  terms  forecast  and  predict  (and  their  derivatives)  are  commonly  used  inter¬ 
changeably  in  meteorology  to  imply  a  projection  fi:t>m  an  initial  time  to  a  later  valid  time. 
This  is  a  prognostic  process.  It  may  employ  d3mamical,  statistical,  subjective,  or  other  tech¬ 
niques.  Unfortunately,  both  terms  are  also  commonly  used  interchangeably  in  statistics 
with  an  entirely  different  connotation.  In  statistics  they  imply  existence  of  an  assodalion 
between  a  dependent  variable  and  a  set  of  independent  (or  explanatoiy)  variables.  This 
relationship  may  be  prognostic,  or  it  may  be  diagnostic  (i.e.,  without  a  time  projection). 

In  statistical  models  such  as  regression,  the  dependent  variable  is  called  the  re¬ 
sponse  while  the  independent  variables  are  called  predictors.  The  terms  dependent  and 
independent  are  unfortunate.  While  meaningful  in  the  vocabulary  of  mathematical  func¬ 
tions,  they  lead  to  confiision  in  this  context  since  the  independent  variables  are  usually 
not  statistically  independent.  Thus,  we  choose  to  use  the  terms  predictor  and  response 
rather  than  dependent  and  independent  variables. 
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Tb  minimize  likely  confusion  resulting  from  a  dual  use  of  the  term  predict,  we  adopt 
the  following  conventions.  The  term  predict  (and  its  derivatives)  will  always  be  used  to 
imply  application  of  a  statistical  model.  As  noted  above,  the  statistical  association  may  be 
either  prognostic  or  diagnostic.  The  term  forecast  (and  its  derivatives)  will  always  be  used 
to  imply  a  prognostic  process,  whether  by  statistical,  dynamical,  or  other  means.  Vi^th 
these  conventions,  a  forecast  may  or  may  not  be  a  prediction,  and  vice  versa. 

1.4.2  Three  Basic  Approaches 

The  three  basic  approaches  to  the  statistical  forecasting  of  weather  parameters  such 
as  cloud  amount  for  a  specified  location,  at  a  specified  valid  time  are:  the  classical  approach, 
the  perfect  prog  approach,  and  the  model  output  statistics  approach.  In  addition,  various 
combinations  and  extensions  are  possible.  The  three  approaches  are  described  below. 

Classical  Approach — ^Predictors  are  specified  at  the  initial  time.  Cloud  amount 
is  forecast  for  the  valid  time  using  a  statistical  method.  Relevant  techniques  include  cor¬ 
relation,  regression,  conditional  probability,  and  autoregressive  moving  average  (ARMA) 
time  series  models.  A  major  shortcoming  of  the  classical  approach  is  that  it  ignores  the 
benefits  of  available  dynamical  forecast  models. 

Perfect  Prog  Approach — ^Predictors  are  specified  at  the  valid  time  by  a  prognos¬ 
tic  process,  usually  numerical.  Cloud  amount  is  then  diagnosed  by  a  statistical  technique. 
Parameters  of  the  diagnostic  model  (e.g.,  the  regression  coefficients)  are  derived  from  ob¬ 
servations,  not  from  forecast  model  outputs.  Since  the  statistical  model  is  unaware  of  the 
forecast  model  biases,  it  cannot  correct  for  them.  From  one  point  of  view,  however,  this 
shortcoming  is  an  advantage:  the  diagnostic  model  will  not  require  revision  when  the  fore¬ 
cast  model  undergoes  a  significant  upgrade. 

Model  Output  Statistics  Approach — Predictors  are  specified  at  the  valid  ume 
by  a  numerical  forecast  model.  Cloud  amount  is  diagnosed  by  a  statistical  model.  Parame¬ 
ters  of  the  diagnostic  model  are  based  on  both  observations  of  cloud  amount  and  on  numer¬ 
ical  model  forecasts  of  the  predictors.  Since  the  statistical  model  coefficients  are  derived 
from  forecast  model  outputs,  the  MOS  forecast  can  correct  inherent  numerical  model 
biases.  The  model  output «  -.tis^ics  approach  usually  outperforms  both  the  classical  ap¬ 
proach  and  perfect  prog  approach  because  it  uses  more  information  than  its  competitors. 
This  assumes,  of  course,  that  an  adequate  development  database  is  available. 
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A  key  point  to  keep  in  mind  concerning  the  MOS  approach  is  that  although  it  can 
account  for  systematic  errors  in  the  driving  model,  it  cannot  account  for  random  model 
errors.  In  addition,  any  significant  changes  to  the  underlying  forecast  model  will  require 
that  a  new  set  of  model  equations  be  developed  based  on  at  least  one  or  two  years  of  model 
output.  In  the  end,  of  course,  the  quality  of  the  MOS  model  is  strongly  dependent  on  both 
the  consistent  performance  and  accuracy  of  the  driving  forecast  model. 

Various  extensions  to  and  combinations  of  the  three  basic  approaches  to  statistical 
cloud  forecasting  are  possible.  We  used  a  combination  of  the  classical  and  model  output 
statistics  approaches  in  this  effort.  That  is,  we  included  both  observations  at  the  initial 
time  and  forecast  variables  at  the  valid  time  in  the  pool  of  potential  predictors. 


1.4^  Recent  Work  Using  the  Model  Output  Statistics  Approach 

The  body  of  literature  on  model  output  statistics  is  very  large.  Here  we  outline  some 
of  the  developments  in  MOS  cloud  amount  forecasting. 

Model  output  statistics  methods  for  cloud  amount  forecasting  date  back  to  the  early 
1970s  and  the  pioneering  work  by  Glahn  at  the  National  Weather  Service’s  Technique  De¬ 
velopment  Laboratory  (Ref  13).  With  Lowry,  he  developed  multiple  linear  regression 
equations  for  coded  cloud  amount  (0=clear,  Ispartial  obscuration,  2=thin  scattered,  . . ., 
7=overcast,  8=obscured)  at  four  stations  for  four  times  of  day.  They  reported  success  using 
observed  cloud  amount,  and  forecasts  of  mean  relative  humidity,  saturation  deficit,  and 
1000  mb  zonal  and  meridional  wind  components.  Interestingly,  in  all  cases,  the  forecast 
time  that  was  selected  for  humidity  variables  by  their  screening  procedure  was  later  than 
the  MOS  forecast  valid  time.  By  using  the  MOS  methodology,  they  were  able  to  detect  a 
phase  error  in  the  model  humidity  field. 

A  defect  in  their  scheme  was  use  of  a  coded  cloud  amount  rather  than  a  categorical 
cloud  amount.  By  1974,  Glahn  realized  this  error  (Refs.  14  and  15)  and  switched  to  a 
REEP  (Regression  Estimation  of  Event  Probability)  forecast  of  categorical  cloud  amoimt 
(see  Section  3. 1  for  a  brief  description  of  the  REEP  technique).  He  also  dropped  the  surface 
observation  as  a  predictor  and  reported  on  a  problem  with  choosing  the  category  with 
highest  predicted  probability,  namely,  biased  forecasts.  By  1976  (Ref  7),  Carter  and 
Glahn  reported  a  solution  to  the  bias  problem:  inflation  of  variance  followed  by  multiplica¬ 
tion  by  the  “minimum  bias  matrix”  adjustment  factor. 
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In  a  1978  AFGL  technical  report  (Ref.  16),  Glahn  reported  on  a  design  of  an  auto¬ 
mated  MOS  system  for  the  Air  Weather  Service.  In  1985,  Perrone  and  Miller  (Ref.  21)  re¬ 
ported  that  the  generalized  exponential  Markov  (GEM)  model  compared  favorably  with 
MOS  for  cloud  amount  forecasts.  In  1988,  Brunet  et  al.  (Ref.  4)  reported  that  for  short- 
range  cloud  amount  forecasts  in  Canada,  perfect  prog  was  superior  to  model  output  statis¬ 
tics  when  the  surface  observations  were  used  as  predictors.  In  1989,  the  opposite 
conclusion  was  reached  by  Carter  et  al.  (Ref.  8)  for  U.S.  forecasts. 

Finally,  in  1990,  Jacks  et  al.  (Ref.  17)  reported  on  a  new  MOS-based  cloud  amount 
forecasting  system  based  on  output  from  the  NWS  Nested  Grid  Model.  Again,  categorical 
cloud  amount  was  forecast,  but  a  new  threshold  technique  was  used  to  remove  bias.  Sig¬ 
nificant  predictors  included  forecasts  of  mean  relative  humidity,  relative  humidity  at  var¬ 
ious  levels,  vertical  velocity,  and  850  mb  moisture  convergence.  In  addition,  Jacks  et  al. 
used  a  variety  of  “binaiy”  predictors,  i.e.,  predictors  whose  value  is  one  if  a  certain  criteri¬ 
on  is  met  and  zero  otherwise.  An  example  of  a  binary  predictor  is  the  relative  humidity 
at  the  surface  tested  against  various  percentage  thresholds. 
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DATA  SOURCES 


Forecasts  from  the  Air  Force  Global  Spectral  Model  and  Real-Ume  Nephanalysis 
data  were  selected  by  the  Government  for  use  during  model  development  in  this  effort. 
The  RTNEPH  data  were  also  used  as  "ground  truth”  because  they  represent  a  unique,  uni¬ 
formly-derived  representation  of  global  cloud  conditions.  In  addition  to  those  two  meteo¬ 
rological  data  sources,  we  used  a  terrain  database  compiled  at  Phillips  Laboratory  on  the 
eighth-mesh  grid  that  contains  terrain  type  (land,  water,  coast)  and  elevation  in  meters 
at  each  grid  point. 

A  large  part  of  oiir  MOS  effort  was  spent  obtaining,  preparing,  and  validating  our 
year-long  developmental  database.  GSM  data  were  first  generated  and  then  remapped  to 
the  Rn^PH  grid  system  by  government  personnel  at  Phillips  Laboratory.  RTNEPH  data 
were  obtained  directly  through  the  Air  Force  Environmental  Tbdinical  Applications  Cen¬ 
ter  (USAFETAC).  Descriptions  of  the  RTNEPH  and  GSM  data  sources  are  provided  in  the 
following  sub-sections.  Along  with  those  descriptions,  we  also  document  the  data  process¬ 
ing  and  validation  tasks  that  were  an  early  part  of  the  effort.  All  the  figures  for  this  section 
are  included  at  the  end  of  the  section,  starting  on  page  19. 

2.1  RTNEPH  GRIDDED  CLOUD  COVER 

2.1.1  Data  Description 

The  RTNEPH  data  archive  has  been  maintained  by  the  U.S.  Air  Force  meteorolo¬ 
gists  since  1983  when  the  RTNEPH  format  replaced  that  of  the  previous  3DNEPH  ar¬ 
chive.  Formats  and  procedures  of  RTNEPH  are  summarized  thoroughly  in  Ref.  23,  but 
documentation  cannot  capture  the  extreme  complexity  of  the  procedures  and  the  d^ree 
of  human  interaction  required  to  produce  this  massive  data  set  operationally  on  a  3-hour 
cycle  without  interruption.  The  RTNEPH  archive  is  produced  by  Air  Force  meteorologists 
at  GWC  and  at  USAFETAC’s  Operating  Location  A  at  Asheville,  NC.  A  complete  global 
RTNEPH  analysis  is  provided  every  three  hours  in  the  northern  hemisphere  and  every 
six  hours  in  the  southern  hemisphere  using  data  from  two  sources:  the  DMSP  satellites 
and  the  Automated  Weather  Network  (AWN). 

AWN  data  include  surface,  upper  air,  and  aircraft  pilot  report  information.  These 
point  data  are  of  great  value  when  and  where  available.  The  fractional  cloud  cover  esti¬ 
mates  by  trained  surface  observers  can  be  especially  valuable  in  establishing  the  output 
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cloud  cover  information  described  below.  Unfortunately,  in  the  U.S.  even  these  data  are 
less  valuable  than  they  might  be  because  Airways  codes  (clear,  scattered,  broken,  and 
overcast)  are  used  by  surface  observers  instead  of  the  "eighths”  categories  used  in  most 
locations  world-wide.  Since  surface  observers  and  reporting  aircraft  are  not  available 
around  the  clock  or  around  the  world  to  support  RTNEPH  requirements,  it  is  on  the  satel¬ 
lite  data  that  most  of  the  data  archive  depends. 

The  Operational  Linescan  Scanners  (OLS)  on  two  DMSP  satellites  (F-8  and  F-9) 
currently  provide  reflectance  and  radiance  data  with  approximately  0.6  km  visible  and 
thermal  IR  resolution  as  they  trace  101  minute,  sun-S3mchronous  orbits  around  the  earth. 
A  third  satellite  (F-10)  is  in  a  non-circular  orbit  and  is  used  only  occasionally  at  GWC. 
Most  points  on  the  earth  are  scanned  four  or  more  times  per  day  by  the  two  OLS  instru¬ 
ments.  Using  thresholding  techniques,  cloud/no-cloud  analyses  are  derived  from  the 
DMSP  data  and  the  fraction  of  cloudy  pixels  within  the  area  represented  by  each  grid 
point  becomes  the  satellite-based  estimate  of  cloud-cover  fraction.  To  the  extent  possible, 
layered  cloud  information  is  extracted  from  the  IR  signals. 

Figure  2  shows  the  importance  of  the  time  flag  reported  in  the  RTNEPH  data.  Due 
to  the  geometry  of  the  satellite  paths,  certain  areas  of  the  globe  are  updated  more  often 
than  others;  the  polar  regions  receive  much  more  coverage  than  the  equatorial  regions. 
The  upper  curve  in  Figure  2  represents  the  time  between  consecutive  DMSP  observations 
as  a  function  of  latitude  along  a  fixed  meridian  (30”  E  longitude)  by  the  ciirrent  constella¬ 
tion.  As  we  will  see  in  Section  2.1.3,  the  age  of  the  data  played  an  important  role  in  the 
selection  of  model  regions  and  times  for  which  to  perform  the  model  development  as  well 
as  in  the  predictor  selection  itself 

Probably  the  most  critical  step  in  the  RTNEPH  procedure  is  that  of  merging  the  two 
data  sources  to  produce  the  best  possible  assessment  at  each  grid  point.  For  many  grid 
points  the  satellite-derived  values  are  the  only  data  available.  In  some  cases,  human  in¬ 
tervention,  called  "bogusing,”  can  occur  during  RTNEPH  generation  especially  at  earth 
locations  of  particular  interest  or  known  to  exhibit  particular  systematic  problems  (e.g., 
coastal  areas  with  persistent  low  clouds,  snow-covered  ground,  etc.). 

Each  eighth-mesh  RTNEPH  analysis  produces  up  to  4096  grid  points  within  each  of 
120  RTNEPH  boxes  (60  in  each  of  the  northern  and  southern  hemispheres).  Grid  point  spac¬ 
ing  varies  as  a  function  of  latitude  from  51  km  near  the  poles  to  26  km  near  the  equator. 
Figure  3  shows  the  60  RTNEPH  boxes  that  make  up  the  northern  hemisphere.  In  the  figure, 
we  highlighted  the  20  boxes  of  interest  in  this  effort  as  selected  by  government  personnel. 
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Data  items  at  each  grid  point  include  both  quantitative  information  regarding 
cloud  cover  and  qualitative  information  on  weather  conditions  and  data  quality  (see 
Table  1).  We  extracted  total  cloud  amount  information,  the  time  flag,  and  eight  diagnostic 
flags  for  our  MOS  development  as  described  below. 

2.1.2  Data  Preparation 

All  data  processing  was  performed  on  the  computers  at  Phillips  Laboratory.  Com¬ 
plete  RTNEPH  analysis  fields  were  received  from  USAFETAC  on  9-track  tape  in  box-time 
order.  Each  tape  holds  four  months  of  data  (reported  every  three  hours:  0,  3, 6, 9, 12, 15 
18,  21Z)  for  one  box.  The  complete  data  set  fits  on  sixty  9-track  tapes  (20  boxes,  for  12 
months).  By  selecting  only  certain  fields  from  the  total  data  set,  we  were  able  to  reduce 
the  number  of  tapes  required  for  this  feasibility  study  to  seven. 

Table  1  RTNEPH  Gridpoint  Information 

Cloud  Data 

Total  Cloud  Amount  (in  percent  of  grid  area  covered) 

For  up  to  four  cloud  layers: 

Cloud  Amount  (rounded  up  to  next  higher  5  percent) 

Cloud  Type 
Cloud  Base  Height 
Cloud  Top  Height 

Weather  Data 

Present  Weather 

Surface  Visibility 

Data  Quality  Information 

Time  Flag  for  Data  Age 

Source  Flag  for  Layer  Data 

Diagnostic  Information  ('t)est  reporT  flags,  viewing  conditions,  bogus  flags, . . .) 


The  exact  format  of  the  USAFETAC  tapes  is  described  in  the  RTNEPH  User’s 
Handbook  (Ref  23).  Briefly,  the  RTNEPH  database  consists  of  a  number  of  variables  per¬ 
taining  to  cloud  cover  (layer  amount,  type,  and  base  and  top  heights  as  well  as  total  cloud 
amount),  qualitative  information  on  weather  conditions,  and  data  quality  indicators.  For 
this  efibrt  we  were  only  concerned  with  total  cloud  amount,  to  be  used  both  as  a  predictor 
and  as  ground  truth.  We  did  not  take  layer  information  into  account  for  a  variety  of  reasons 
going  beyond  our  desire  to  control  costs  and  computational  complexity.  For  example,  it  is 
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known  that  the  KHTNEPH  analysis  model  is  less  accurate  for  middle  layer  cbuds  than  for 
the  highest  and  lowest  layers  which  can  be  observed  by  satellites  and  human  daservers. 
Furthermore,  many  critical  operational  missions  can  be  served  with  total  cloud  amount 
forecasts  alone. 

Along  with  header  information  (including  box  number,  year,  month,  day,  and  hour), 
we  extracted  three  single-byte  variables  for  each  eighth-mesh  grid  point.  They  are: 

•  total  cloud  cover  amount  (0-100%) 

•  time  in  hours  of  the  newest  data  at  the  point  (0-255) 

•  8-bit  diagnostic  value,  where  each  of  the  8  bits  represents  one  diagnostic  flag 
as  follows: 

bit  1  —  point  was  bogused 

bit  2  —  conventional  report  was  included 

bit  3  —  data  from  a  conventional  report  was  spread  to  this  point 

bit  4  —  visual  satellite  data  were  included  in  this  analysis 

bit  5  —  IR  satellite  data  were  included  in  this  analysis 

bit  6  —  low  level  cloud  has  been  persisted  past  normal  data  cutoff  time 

bit  7  —  although  the  visible  satellite  observed  cloudy,  point  was  marked 
clear  for  lack  of  supporting  data 

bit  8  —  fog/haze  superseded  by  other  weather  elements  in  present  weather. 

As  we  will  discuss  later  in  Section  3,  we  used  the  total  cloud  amount  directly  as  a 
predictor  and  as  ground  truth.  The  data  flags  were  used  to  eliminate  old  and/or  question¬ 
able  data  points  from  the  regression  analysis. 

2.1^  Data  Validation 

RTNEPH  data  validation  occurred  during  tape  processing  at  Phillips  Laboratory 
and  later  on  TASC  computers.  Our  goals  during  this  task  went  beyond  simply  validating 
the  data,  to  include  data  analysis.  In  validating  the  RTNEPH  data  we  were  able  to  detect 
patterns  in  the  data  and  consider  the  effects  of  data  timeliness,  cloud  cover  distributions, 
and  artifacts  on  the  regression  model  development.  Data  validation  and  analysis  con¬ 
sisted  of  the  following  components: 

•  check  headers  for  correctness 

•  verify  that  data  fields  in  adjacent  boxes  are  continuous  across  box  boundaries 

•  compare  total  cloud  amount  fields  with  daily  weather  maps 
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•  visually  analyze  concurrent  cloud  a>ver,  time,  and  diagnostic  flags 

•  generate  histograms  of  cloud  cover,  time,  and  diagnostic  flags. 

The  last  three  components  deserve  further  description. 

Comparison  with  weather  data — We  compared  cloud  cover  fields  in  box  45, 
which  includes  the  U.S.  east  coast,  with  daily  weather  maps  produced  by  the  NOAA  Cli¬ 
mate  Analysis  Center  for  the  same  region  over  a  three-month  period  (February-April 
1989).  Figure  4  shows  an  example  in  which  a  strong  front  is  present  in  the  daily  weather 
map  over  the  east  coast  stretching  from  '^rginia  to  Nova  Scotia.  The  clouds  associated 
with  that  front  are  clearly  visible  in  the  corresponding  RTNEPH  data  field.  Over  the 
three-month  period  we  were  able  to  qualitatively  verify  the  RTNEPH  analysis  fields  on 
days  that  had  strong  cloud  signals. 

Visual  analysis — In  this  task  we  simultaneously  displayed  (for  a  selected  box/ 
time)  gridded  cloud  cover  amount,  data  age,  and  eight  diagnostic  flags.  By  viewing  all  ten 
data  fields  at  one  time,  we  were  able  not  only  to  see  artifacts  in  the  cloud  amount  fields, 
but  also  to  use  the  ancillaiy  diagnostic  data  to  determine  the  source  of  the  artifacts.  For 
example,  linear  artifacts  are  commonly  found  in  RTNEPH  cloud  amount  fields.  By  ana¬ 
lyzing  the  time  and  diagnostic  flags,  it  is  clear  that  those  linear  artifacts  frequently  sepa¬ 
rate  areas  of  distinct  "age”  and  form  the  boundaries  of  satellite  swaths  through  the  region. 
Figure  5  shows  an  example  of  just  this  tjrpe  of  artifact. 

The  upper-left  image  in  Figura  5  shows  fi’actional  cloud  cover  (black  represents 
clear,  white  is  used  for  overcast,  and  a  gray  scale  is  used  to  represent  partially  cloudy  grid 
points).  The  next  image  to  the  right  shows  the  relative  age  of  the  data  used  in  the 
RTNEPH  analysis.  Here  the  gray  scale  changes  from  dark  for  newer  data  to  light  for  older 
data.  The  next  eight  binary  images  are  the  eight  diagnostic  flags.  In  those  images  white 
represents  0  (the  flag  is  not  set)  and  black  represents  1  (the  flag  is  set).  The  8  diagnostic 
flags  are  described  briefly  in  Section  2.1.2. 

Other  common  features  foimd  in  the  RTNEPH  cloud  amoiuit  fields  are  the  star-  and 
cross-shaped  artifacts  that  are  a  sign  of  "spread”  data.  That  is,  doud  cover  at  one  grid  point 
is  "spread”  into  atyacent  grid  points.  Again,  by  viewing  the  diagnostic  flags,  data  age,  and 
doud  cover  fields  simultaneously,  we  were  able  to  determine  that  this  type  of  artifact  is  a 
characteristic  of  the  RTNEPH  analysis  algorithms,  not  a  data  processing  error.  Figure  6 
shows  a  case  where  spread  data  are  clearly  visible.  The  format  of  the  ten  images  is  the 
same  as  described  above.  We  dedded  to  remove  all  spread  data  points  from  our  model  de¬ 
velopment  data  set. 
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The  visual  analysis  task  was  important  not  only  to  verify  our  data  processing,  but  to 
understand  the  strengths  and  weaknesses  of  the  RTNEPH  data  fields  and  the  effects  of  both 
on  regression  model  development.  By  viewing  a  large  number  of  data  fields  we  were  able  to 
qualitatively  compare  data  from  different  regions,  seasons,  and  time-of-day.  We  looked  for 
the  number  and  type  of  artifacts  present,  the  relative  percentage  of  surface  observations 
and  satellite  data  used,  the  age  of  the  data  and  how  it  changed  with  time  of  day,  the  fi^uen- 
cy  of  bogused  data  points  and  any  missing  data.  From  this  analysis,  we  were  able  to  esti¬ 
mate  the  potential  limitations  of  the  RTNEPH  data  both  as  a  predictor  and  ground-truth 
data  source,  especially  with  respect  to  the  age  of  the  data  in  certain  regions. 

Histogram  generation — ^This  was  the  most  quantitative  and  telling  of  the  data  val¬ 
idation  and  analysis  tasks.  Following  directly  fi^om  the  qualitative  visual  analysis  described 
above,  we  generated  histograms  for  the  cloud  cover,  time  flag,  and  two  of  the  diagnostic  flags 
to  better  imderstand  these  data  and  their  potential  effects  on  model  development.  Figures  7 
through  15  include  histograms  of  fiactional  cloud  cover  and  the  data  time  flag  for  several 
boxes,  seasons,  and  times  that  are  representative  of  the  many  cases  we  looked  at.  These  fig¬ 
ures  are  discussed  below  (refer  to  Figure  3  for  the  geographic  locations  of  each  box). 

We  begin  with  a  discussion  about  the  time  flag  and  its  significance.  A  key  element 
of  short-range  (temporal)  forecasting  is  up-to-date  initialization  data.  A  key  element  of 
forecast  validation  and  regression  modeling  is  up-to-date  validation  data.  Because  we 
used  the  RTNEPH  cloud  cover  fields  in  both  roles,  its  temporal  characteristics  were  espe¬ 
cially  important  for  model  development. 

In  many  regions  of  the  globe,  especially  in  the  tropics  and  over  water,  the  RTNEPH 
is  almost  solely  driven  by  satellite  observations  at  sparse  intervals,  thus  there  are  large 
gaps  in  the  data.  These  data  gaps  limited  the  number  of  regions  and  times  for  which  we 
could  develop  meaningful  models  and  from  which  we  could  draw  meaningful  conclusions 
regarding  forecast  quality  based  on  comparisons  with  persistence. 

Figure  7  shows  the  fraction  of  cloud  data  in  the  RTNEPH  that  is  less  than  three 
hours  old  for  box  34  (mid  Pacific).  This  plot  is  representative  of  the  problems  in  ocean  re¬ 
gions:  very  few,  if  any,  automated  weather  sites  and  scanty  satellite  coverage,  resulting 
in  large  time  gaps  between  cloud  analysis  updates.  A  plot  of  the  percentage  of  new  (i.e., 
less  than  three  hours  old)  data  in  box  61  (including  Colombia  and  Venezuela),  shown  in 
Figure  8,  shows  similar  characteristics.  Figure  9  shows  a  histogram  of  the  data  age  for  box 
45  (east  coast  U.S.  and  western  Atlantic)  at  21Z.  In  this  example,  less  than  20%  of  the  data 
are  less  than  7  hours  old. 
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On  the  other  hand,  box  30,  which  includes  eastern  Europe  and  western  Russia,  has 
a  large  automated  weather  network  and  more  frequent  satellite  coverage.  Thus,  a  large  per¬ 
centage  of  the  grid  points  are  “new,”  as  shown  in  Figime  10.  A  histogram  of  box  44  data  (cen¬ 
tral  U.S.  and  Canada)  shows  similar  results  (see  Figure  11).  Those  cases  with  the  most 
up-to-date  data  will  provide  the  best  evaluation  of  the  potential  value  of  the  MOS  approach. 

We  found  a  wide  variety  in  cloud  cover  means  and  distributions  across  the  20-box 
data  set.  Some  of  that  variety  can  be  seen  in  Figures  12  through  15.  For  example.  Figures  12 
and  13  (containing  data  for  box  61  afternoon  and  evening,  respectively)  show  the  diurnal 
differences  in  cloud  cover  characteristics  of  the  tropics.  Figure  12  has  almost  15%  fewer  grid 
points  in  the  two  end  categories,  clear  and  overcast,  than  Figure  13.  Thus,  a  greater  per¬ 
centage  of  the  grid  points  in  the  afternoon  are  partly  cloudy  which  is  characteristic  of  pop¬ 
corn  cumulus  which  depend  strongly  on  diurnal  heating. 

Cloud  patterns  in  the  winter  at  midlatitudes  (Figures  14  and  15)  show  a  different  diur¬ 
nal  signature  whose  dominant  feature  is  the  pronoimced  clearing  at  ni^t.  A  unique  artifact 
present  only  in  the  those  boxes  containing  regions  of  the  U.S.  is  the  peaks  in  doud  cover  at  25 
and  75%.  These  two  values  correspond  to  the  U.S.  Airwa3rs  codes  for  scattered  and  broken,  re¬ 
spectively,  and  can  be  seen  dearly  in  these  two  figures  corresponding  to  doud  cover  in  box  44. 

By  analyzing  these  histograms  and  many  others  we  were  able  to  quantitatively 
verify  the  RTNEPH  data.  In  addition,  this  analysis  led  us  to  concentrate  our  model  devel¬ 
opment  efforts  on  those  regions  and  times  for  which  we  would  best  be  able  to  evaluate  the 
MOS  approach  with  respect  to  persistence  for  3-,  6-,  and  9-hour  forecasts.  This  precluded 
those  regions  and  times  for  which  little  up-to-date  cloud  cover  data  were  available,  be¬ 
cause  in  those  cases,  persistence  dominated.  As  an  extreme  example,  if  all  the  grid  points 
for  a  given  valid  time  in  a  development  data  set  had  cloud  amounts  that  had  not  been  up¬ 
dated  in  over  9  hours,  a  comparison  of  the  valid-time  RTNEPH  cloud  cover  to  the  initial¬ 
time  (persisted)  cloud  cover  would  result  in  100%  skill  for  3-,  6-,  and  9-hour  forecasts. 
Similarly,  the  MOS-produced  forecasts  would  perform  100%  because  they  also  use  per¬ 
sisted  cloud  cover  in  the  predictor  pool.  Evaluating  these  forecasts  would  hardly  provide 
a  worthwhile  or  meaningful  assessment  of  the  MOS  approach. 

On  the  other  hand,  though  evaluation  requires  only  the  most  current  data,  operation¬ 
al  forecast  models  must  use  whatever  data  are  available  at  initialization  time.  Keeping  both 
of  these  requirements  in  mind,  we  included  only  those  data  points  whose  valid-time 
RTNEPH  age  was  less  than  three  hours  in  model  development  and  validation,  but  we  had 
no  cutoff  age  for  the  persisted  RTNEPH  cloud  cover  used  as  a  predictor  in  the  MOS  model. 
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2.2  GLOBAL  SPECTRAL  MODEL  OUTPUT 


2J2.1  Data  Description 

The  Air  Force  Phillips  Laboratoi7  maintains  and  runs  a  copy  of  the  twelve-layer 
Global  Spectral  Model  that  was  the  primary  source  of  model  data  for  this  effort.  The  GSM 
(Ref.  24)  is  based  on  integration  of  equations  representing  five  prognostic  variables; 

absolute  vorticity 

divergence 

temperatiure 

surface  pressure  (logarithm) 
specific  humidity. 

From  these,  four  additional  diagnostic  variables  are  computed  routinely:  geopotential  height, 
vertical  velocity,  and  the  u  (zonal)  and  v  (meridional)  wind  components. 

The  vertical  coordinate  used  within  the  model  is  the  normalized  pressure  or  "sigma” 
coordinate.  Typically,  model  variables  are  interpolated  to  constant  pressure  mandatory  levels 
for  output  as  shown  in  Table  2.  However,  in  this  effort  we  were  provided  with  wind,  tempera¬ 
ture,  and  humidity  variables  for  sigma  layers  as  used  in  the  model,  thus  avoiding  errors 
introduced  by  interpolating  from  sigma  to  mandatory  levels.  The  thirteen  sigma  levels  that 
define  the  twelve  sigma  layers  used  in  the  GSM  are  shown  in  Figure  16. 

Moisture  is  represented  in  the  GSM  by  the  mixing  ratio,  or  equivalently,  specific 
humidity  (one  of  the  prognostic  model  variables).  Evaporation  firom  the  oceans  is  repre¬ 
sented  in  the  model  including  the  impact  of  the  surface  wind  speed  on  the  evaporation 
rate.  Evaporation  from  land,  however,  is  not  represented  and  is  a  known  shortcoming  of 
the  present  model. 

Spatial  resolution  of  the  GSM  is  limited  by  the  40-wave  truncation  of  the  spherical 
harmonics  used  to  represent  the  meteorological  fields  at  each  pressure  level.  The  40- wave 
model  can  support  resolutions  down  to  250  km  near  the  equator,  but  it  is  often  evaluated 
on  a  2.5**  by  2.5°  grid  (about  280  km  x  280  km  at  the  equator).  For  this  effort  we  obtained 
GSM  outputs  on  the  half-mesh  grid  (approximately  100  km  grid  spacing  near  the  equator) 
that  were  produced  by  interpolating  from  the  2.5°  by  2.5°  grid  to  the  half-mesh  grid.  That 
interpolation  was  performed  by  personnel  at  Phillips  Laboratory. 
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Table  2  Air  Force  GSM  Forecast  Fields 


FIELD 

SURFACE 

1 

1 

O 

1000  mb 

850  mb 

700  mb 

500  mb 

400  mb 

300  mb 

250  mb 

200  mb 

150  mb 

100  mb 

1000-650  mb 

Sea-level  Pressure 

X 

Heights  (D-Value) 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

Temperature 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

Temperature  Advection 

X 

U-component  Wind 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

V-component  Wind 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

Vertical  Velocity 

X 

X 

X 

X 

X 

X 

X 

X 

X 

Relative  Humidity 

X 

X 

X 

X 

X 

Specific  Humidity 

X 

X 

X 

X 

X 

X 

Dewpoint  Depression 

X 

X 

X 

X 

X 

Vorticity 

X 

X 

X 

Vorticity  Advection 

X 

Stream  Function 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

Tropopause  Pressure 

Tropopause  Temperature 

Tropopause  Height 

The  GSM  generates  forecasts  every  15  minutes  out  to  24  hours  and  is  initialized 
twice  a  day  at  OZ  and  12Z.  For  this  effort  we  obtained  forecasts  every  3  hours  to  match  the 
temporal  resolution  of  the  RTNEPH  data.  We  obtained  both  the  OZ-  and  the  12Z-initial- 
ized  model  runs  eveiy  day  of  1989  for  the  twenty  boxes  shown  in  Figure  3.  For  this  study 
GSM  model  runs  were  initialized  with  historical  HIRAS  data  and  thus  were  identical  to 
operational  model  runs  at  Global  Weather  Central  for  that  same  period. 

2.2JS  Data  Preparation 

The  Global  Spectral  Model  runs  were  performed  on  the  Cray2  at  the  Weapons  Lab¬ 
oratory  in  Albuquerque,  NM  by  Phillips  Laboratory  personnel.  In  a  post-processing  step, 
spectral  outputs  of  the  GSM  were  converted  to  latitude-longitude  coordinates,  vorticity 
and  divergence  values  were  translated  to  u  and  v  wind  components,  and  specific  humidity 
was  transformed  to  relative  humidity.  Finally,  a  nine-point  quadratic  interpolation  was 
used  to  interpolate  from  latitude-longitude  coordinates  to  the  half-mesh  RTNEPH  grid. 
These  data  were  then  passed  on  to  TASC  for  further  processing. 
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As  part  of  our  data  processing  task,  we  sorted  the  processed  GSM  outputs  from  syn¬ 
optic  order  to  box/time  order  using  the  computers  at  Phillips  Laboratory.  We  used  box/ 
time  arrangement  throughout  this  effort  because  it  was  the  most  convenient  for 
developing  regional  models.  During  the  sorting  process  we  also  scaled  the  32-bit  data  val¬ 
ues  at  each  grid  point  to  16-bit  integer  values  to  save  space  and  speed  up  processing,  and 
yet  retain  the  desired  accuracy. 

The  resulting  data  files  contained  a  header  for  each  forecast  time  that  included:  box 
number,  year,  month,  day,  model  initialization  time  (OZ  or  12Z),  forecast  time  (3-,  6-, . . ., 
24-hour  forecasts),  and  sigma  layer  heights.  The  remainder  of  the  data  files  contained  the 
following  gridded  variables  on  the  half-mesh  grid  at  up  to  ten  vertical  levels:  u  wind  (10 
sigma  levels),  v  wind  (10  sigma  levels),  temperature  (10  sigma  levels),  relative  humidity 
(lowest  six  sigma  levels),  surface  pressure,  and  surface  height. 

Because  our  model  development  data  set  contains  model  runs  initialized  both  at  OZ 
and  12Z,  and  because  each  run  produces  forecasts  out  to  24  hoiurs,  we  have  two  forecasts 
available  at  any  given  vaUd  time.  For  example,  for  a  3Z  valid  time,  we  have  the  3-hour 
forecast  from  the  OZ  model  run  of  the  same  day  and  the  15-hour  forecast  from  the  12Z  mod¬ 
el  run  of  the  previous  day.  However,  both  of  these  forecasts  are  not  always  available  due 
to  operational  schedules.  In  otu*  example  of  a  3Z  valid  time,  output  from  the  model  run 
initialized  at  OZ  would  not  be  immediately  available  at  OZ  in  order  to  make  a  3-hoiu'  fore¬ 
cast  valid  at  3Z.  Similarly,  there  would  be  no  OZ-initialized  data  available  for  the  6-  or 
9-hour  forecasts  valid  at  3Z.  Figures  17  through  20  show  time  lines  for  each  of  four  fore¬ 
cast  valid  times  (0, 3, 6,  and  9Z)  and  the  GSM  data  available  for  the  3-,  6-,  and  9-hoiu’  fore¬ 
casts  at  each  time.  Similar  timelines  can  be  made  for  forecast  times  12Z  through  21Z  in 
which  the  OZ  and  12Z  GSM  data  sources  are  switched. 

For  any  given  box,  season,  and  valid  time,  a  single  model  development  data  file  was 
built  that  included  the  GSM  data  for  the  valid  time,  the  corresponding  RTNEPH  cloud 
amount  and  diagnostic  fields,  along  with  RTNEPH  cloud  cover  from  the  forecast  initial¬ 
ization  time  back  to  24  hours  before  the  valid  time  to  use  for  persistence. 

2.2.3  Data  Validation 

Data  validation  tasks  for  the  Global  Spectral  Model  data  fields  were  similar  to 
those  we  performed  for  the  RTNEPH  data.  Briefly  they  were: 

•  check  headers  for  correctness 

•  verify  that  data  fields  in  adjacent  boxes  are  continuous  across  box  boundaries 
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•  visuaUy  analyze  data  fields 

•  generate  hisU^prams  of  model  variables. 

The  last  two  tasks  are  described  in  more  detail  below. 

Visual  analysis — ^In  this  task  we  displayed  the  u,  v,  temperature,  relative  humid¬ 
ity,  siuface  pressure,  and  surface  height  fields  for  many  selected  boxes  and  times.  We  dis¬ 
played  all  layer  data  at  one  time  to  verify  (X>ntinuity  in  the  vertical  dimension.  All  data 
were  displayed  as  gray-scale  images.  We  analyzed  the  images  for  potential  artifacts  and 
overall  structure.  We  did  not  find  any  questionable  data  fields. 

Figures  21  through  23  show  examples  of  the  images  analyzed.  The  first  of  those  fig¬ 
ures  shows  the  u  field  for  a  specific  day  and  time  in  box  30.  The  upper-left  image  in  the  se¬ 
quence  is  the  surface-level  u  field  and  each  consecutive  image  shows  the  u  field  at  the  next 
higher  sigma  layer.  We  use  a  gray  scale  to  show  the  relative  magnitude  of  the  field  value, 
where  black  is  low  and  white  is  high.  Figure  22  shows  the  similar  images  for  the  correspond¬ 
ing  temperature  field.  Figure  23  continues  with  the  same  box  and  time  and  includes  the 
relative  humidity  at  the  lower  six  sigma  layers  of  the  atmosphere  and  the  surface  pressure 
and  height  fields.  (Note  how  the  surface  hei^t  field  follows  the  elevation  contours  we  expect 
to  see  for  the  region  in  Europe  corresponding  to  box  30.) 

Histc^ram  generation— As  a  means  to  validate  and  better  understand  the  MOS 
data,  we  generated  histograms  to  show  the  overall  range,  mean,  and  distribution  of  each 
model  variable  in  difierent  locations  and  for  different  seasons  and  times.  This  activity  nat¬ 
urally  led  to  generating  histograms  for  meteorological  quantities  derived  firom  the  MOS 
data  as  well,  such  as  vertical  velocity,  diveigence,  etc.  Example  histograms  of  derived 
quantities  are  shown  in  Section  3.3.  Here,  in  Figures  24  and  25,  we  provide  examples  of 
histograms  of  the  GSM  variables;  u,  v,  temperatiure,  relative  humidity,  and  surface  pres¬ 
sure  for  two  different  box  seasons. 
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Figure  2  Satellite  Weather  Data  Access  Performance — Average  Case  Sample 


19 


NOUtHdUl 
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Figure  4  Coincident  Surface  Weather  and  RTNEPH  Analysis  for  Morning  of 
March  29, 1989  Showing  Similar  Cloud  Cover  Over  the  East  Coast  in 
the  Common  Region  Outlined  in  Black 
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Figiure  7  Histogram  of  the  Percentage  of  Grid  Points  That  are  Less  Than 
3  Hours  Old  as  a  Function  of  Time  of  Day  in  Box  34,  Spring 
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BOX  61.  SUMMER 
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Figure  9  Histogram  Showing  the  Percentage  of  Grid  Points  of 
Various  Ages  in  Box  45  Over  the  Spring  Season  at  2 IZ 
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%  OF  GRIDPOINTS  WITH  "NEW"  DATA 
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BOX  44.  WINTER 
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Figure  11  Histogr2un  of  the  Percentage  of  Grid  Points  That  are  Less  Than 

3  Hours  Old  as  a  Function  of  Time  of  Day  in  Box  44,  Winter 
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BOX  61,  SUMMER.  18Z  (APPROX.  1PM  LOCAL  TIME) 


CLOUD  AMOUMT(%) 

Figure  12  Cloud  Cover  Distribution  for  Box  6 1,  Suzomer 

at  Approximately  1  p.m.  Local  Time 
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BOX  61.  SUMMER.  6Z  (APPROX.  1AM  LOCAL  TIME) 
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Figure  13  Cloud  Cover  Distribution  for  Box  6 1,  Summer 

at  Approximately  1  a.m.  Local  Time 
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Figure  14  Cloud  Cover  Distribution  for  Box  44,  Winter 

at  ^proximately  2  p.m.  Local  Time 
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BOX  44.  WINTER.  9Z  (APPROX.  2AM  LOCAL  TIME) 
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Cloud  Cover  Distribution  for  Box  44.  Winter 
at  Approximately  2  a.m.  Local  Time 
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Figure  17  Schematic  Showing  Whidi  RTNEPH  and  GSM  Data  Are  Available 
for  3-,  6-,  and  9-hour  Forecasts  Valid  at  OZ  (‘x’  Marks  the  Data 
Sources  Used  in  Model  Development) 
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Figure  18  Schematic  Showing  Which  RTNEPH  and  GSM  Data  Are  Available 

for  3-,  6-,  and  9-hour  Forecasts  Valid  at  3Z  (‘x’  Marks  the  Data 
Sources  Used  in  Model  Development) 
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Figure  19  Schematic  Showing  Which  RTNEPH  and  GSM  Data  Are  Available 

for  3-,  6-,  and  9-hour  Forecasts  Valid  at  6Z  (V  Marks  the  Data 
Sources  Used  in  Model  Development) 
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Figure  20  Schematic  Showing  Which  RTNEPH  and  GSM  Data  Are  Available 

for  3-,  6-,  and  9-hour  Forecasts  Valid  at  9Z  (‘x’  Marks  the  Data 
Sources  Used  in  Model  Development) 
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Figure  22  Gray-Scale  Images  Showing  GSM  Forecast  Temperature  at  10  Consecutive  Sigma  Layei 

the  Surface  and  Layer  10  Corresponds  to  Approximately  100  mb)  for  Box  30,  June  1, 15Z 


Figure  23  Gray-Scale  Images  Showing  GSM  Forecast  Relative  Humidity  at  6  Consecutive  Sigma  Layers  (Layer  1  is  the  Surface 

Layer  6  Corresponds  to  Approximately  300  mb),  Surface  Pressure,  and  Surface  Height  Fields  for  Box  30,  June  1, 15Z 
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Figure  24  Histograms  of  Raw  GSM  Variables:  u,v,  Temperature,  Relative  Humidity,  and 

Pressure  at  the  Surface  for  Box  44,  December,  OZ 
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MODEL  DEVELOPMENT  AND  RESULTS 


In  this  section  we  present  the  results  of  our  modeling  effort.  A  general  introduction 
of  various  regression  techniques  and  the  mathematics  behind  developing  linear  regres¬ 
sion  models  can  be  found  in  most  introductory  statistics  texts  (Refs.  10  and  19  are  good 
examples)  and  will  therefore  not  be  included  here.  We  will  introduce  those  concepts 
unique  to  our  model  implementation  and  development  data  set,  but  will  otherwise  as¬ 
sume  that  the  reader  has  a  basic  knowledge  of  linear  regression  for  the  following  discus¬ 
sion  and  the  remainder  of  this  report. 

Throughout  the  model  development  process,  we  employed  a  commercial  statistical 
and  graphical  analysis  software  package  called  STATISTICA~.  It  served  as  a  database 
management  system  for  the  model  development  data  sets  and  provided  a  variety  of  regres¬ 
sion  techniques  and  corresponding  significance  tests,  residual  analysis,  and  graphing  that 
were  critical  to  efficient  modeling  and  analysis  during  this  project.  The  STATISTICAL  soft¬ 
ware  package  provided  modules  to  perform  multiple  linear  regression,  using  least  squares 
estimation,  and  nonlinear  regression  in  which  the  user  can  estimate  any  arbitrary  regres¬ 
sion  function  and  any  arbitrary  loss  function.  It  also  enabled  us  to  save  the  forecast  cloud 
amounts  and  residuals  to  compute  probabilistic  forecasts  off-hne. 

All  of  the  figures  discussed  in  this  section  are  included  at  the  end  of  the  section, 
starting  on  page  96. 

3.1  OVERVIEW  OF  THE  MODELING  PROCESS 

Model  development  consisted  of  the  following  four  major  steps.  These  four  steps  are 
described  at  length  in  subsections  3.2,  3.3, 3.4,  and  3.5,  respectively. 

STEP  1  Build  clusters  based  on  cloud  climatology  —  The  first  step  in  the 
model  development  sequence  was  to  reduce  the  number  of  regression 
models  using  principal  components  analysis.  Regression  models  were 
built  for  each  cluster  rather  than  for  individual  grid  points.  We  per¬ 
formed  cluster  analysis  using  principal  components  for  four  box/season 
combinations  including  one  polar,  one  tropical,  and  two  midlatitude 
boxes.  Variables  in  the  analysis  included:  cloud  cover  mean  and  scale 
distance,  terrain  t)rpe,  elevation,  and  grid  point  location. 
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STEP  2  Build  the  pool  of  potential  predictors — Using  knowledge  of  meteo¬ 
rological  processes  in  general,  cloud  formation  in  particular,  and  some 
trial-and-error,  we  selected  variables  from  GSM  outputs,  RTNEPH 
analyses,  and  the  RTNEPH  terrain  database  to  include  in  the  pool  of  po¬ 
tential  predictors.  In  addition,  we  derived  other  meaningful  quantities 
from  these  outputs  for  inclusion.  We  used  this  common  pool  of  predictors 
for  all  model  development.  We  then  selected  a  subset  of  these  variables 
for  each  box/season/time  using  a  forward  regression  technique  that 
ranked  the  predictors  in  order  of  significance.  Using  only  the  most  signif¬ 
icant  predictors,  we  estimated  regression  coefficients  in  Step  3. 

STEP  3  Develop  regression  models — Many  regression  techniques  are  avail¬ 
able.  These  include  multiple  linear  regression,  nonlinear  regression,  re¬ 
gression  estimation  of  event  probability,  empirical  orthogonal  functions, 
regression  trees,  and  many  others.  We  used  multiple  linear  regression  to 
predict  total  cloud  amount  at  eighth-mesh  resolution  with  the  predictors 
selected  in  the  previous  step.  We  developed  equations  for  selected  cases 
(i.e.,  a  particular  box,  season,  and  time).  We  extended  model  development 
to  predict  probabilistic  forecasts  of  cloud  amount  categories  by  analyzing 
the  distributions  of  the  residuals  as  described  later  in  this  section. 

STEP  4  Evaluate  model  performance — It  is  virtually  always  the  case  that  a 
statistical  forecasting  scheme  performs  better  for  the  development  data 
set  than  for  a  new,  independent  data  set.  This  may  be  due  to  fitting  noise 
in  the  development  data  set  or  to  the  incompleteness  of  the  development 
data  set.  We  used  the  jackknife  technique  to  estimate  the  performance 
of  oiu*  MOS  approach  and  to  determine  confidence  limits  on  that  esti¬ 
mate.  In  addition,  we  also  compared  the  MOS  results  to  those  obtained 
using  simple  persistence.  The  metrics  we  used  in  this  comparison  were 
forecast  sWll  scores  (the  Brier  score,  sharpness,  and  the  20/20  score  dis¬ 
cussed  later  in  this  section).  Lastly,  we  validated  the  probabilistic  fore¬ 
casts  using  Bayesian  analysis  of  the  probability  distributions  obtained 
using  the  model  development  data  set  and  an  independent  set  of  data. 


3.1.1  The  Proposed  Methodology  Using  KEEP 

In  our  original  proposal  (Ref.  26),  we  proposed  developing  MOS  equations  to  fore¬ 
cast  total  cloud  amount  and  probabilities  of  cloud  amount  categories  using  the  regression 
estimation  of  event  probability  methodology  (a  linear  regression  technique)  to  support 
two  different  forecast  formats  required  by  Air  Force  and  Army  decision  makers.  It  turned 
out  that  we  were  able  to  produce  both  forecast  output  types  using  a  more  efficient  and  di¬ 
rect  multiple  linear  regression  technique  which  avoided  some  of  the  pitfalls  of  the  regres¬ 
sion  estimation  of  event  probability  (REEP)  method. 
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The  KEEP  method  (Ref.  22)  requires  a  categorical  response  variable.  The  first  step 
is  to  transform  the  RTNEPH  continuous  cloud  amount  to  a  categorical  variable  by  binning 
cloud  cover  into  N  non-overlapping  categories.  The  single  response  used  in  multiple  linear 
regression  is  then  replaced  with  N  binary  response  variables.  For  each  observation,  the 
response  is  one  if  the  observed  cloud  amount  falls  within  the  corresponding  category  and 
zero  otherwise.  Multiple  models  are  built,  one  for  each  of  the  categories,  all  with  exactly 
the  same  predictor  set  and  exactly  the  same  database.  The  output  of  each  model  is  an  esti¬ 
mate  of  the  probability  that  the  observed  cloud  amount  is  within  the  corresponding  range 
for  a  given  realization  of  the  predictors.  A  REEP  model  is  constructed  to  guarantee  that 
the  sum  of  the  N  probabilities  is  one.  The  predicted  total  cloud  amount  can  be  determined 
from  the  N-vector  of  probabilities  using  a  prescribed  selection  criterion  (e.g.,  the  most 
probable  category,  the  median  category,  etc.). 

REEP  posed  two  problems  in  our  modeling  effort.  First,  it  is  well  known  that  al¬ 
though  the  sum  of  the  resulting  probabilities  is  one,  linear  regression  does  not  constrain 
any  individual  probability  to  lie  between  zero  and  one,  as  a  true  probability  must.  It  is 
possible  to  deal  with  this  problem  by  using  a  nonlinear  response  function  such  as  the  logis¬ 
tic  function  (see  Ref  1),  but  this  increases  the  computational  costs  immensely  and  for  our 
data  set  was  therefore  not  a  feasible  solution.  Alternatively,  it  is  possible  to  scale  the  prob¬ 
abilities,  or  to  simply  set  those  probabilities  less  than  zero  to  zero  and  those  greater  than 
one  to  one.  This  is  rather  an  ad  hoc  and  unsatisfying  solution. 

The  second  problem  is  the  increased  computational  costs  of  developing  the  initial 
predictor  data  set.  There  is  no  reason  to  believe  that  the  same  variables  that  are  strong 
predictors  of  overcast  are  also  strong  predictors  of  clear  or  partly  cloudy.  Therefore  one 
must  determine  the  strongest  predictors  for  each  of  the  cloud  cover  categories  from  the 
overall  pool  of  potential  predictors  individually.  That  is  a  computationally  intensive  pro¬ 
cess  since  the  number  of  potential  predictors  is  large  (see  Section  3.3).  After  selecting  the 
most  significant  predictors  for  each  of  the  categories  one  can  build  a  group  that  is  common 
to  all  of  the  categories  (a  somewhat  subjective  process)  and  perform  the  regression  using 
the  common  predictor  data  set.  The  actual  work  required  then  to  perform  the  N  multiple 
regressions  from  the  same  predictor  set  is  little  more  than  that  to  perform  a  single  regres¬ 
sion  because  the  covariance  matrix  does  not  change. 
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In  Section  3.4.3  we  discuss  the  particulars  of  how  we  implemented  KEEP  and  show 
the  results  of  a  sample  case.  Now  we  continue  with  a  brief  description  of  the  alternate  meth¬ 
odology  we  adopted  to  replace  KEEP,  and  discuss  how  we  generated  both  forecast  output 
types  (total  cloud  amount  and  probabilities  of  cloud  amount  categories)  from  a  single  lin¬ 
ear  regression  model. 

3.1.2  An  Alternate  Methodology  Using  Actual  Residuals 

Our  goal  was  to  find  an  empirical  linear  function  in  which  the  response  variable, 
total  cloud  amount,  is  expressed  in  terms  of  k  predictor  variables  selected  from  GSM  out¬ 
put,  RTNEPH  analyses,  and  other  candidate  predictors.  It  is  important  to  note  that  this 
model  is  linear  in  the  predictors,  but  not  necessarily  in  the  meteorological  variables,  as 
predictors  may  be  transformations  of  meteorological  variables. 

Following  the  notation  in  Ref.  6,  let  xi,  X2, .  . .,  xit  be  k  predictors  and  y  be  the  re¬ 
sponse  variable,  then  the  linear  regression  model  is 

yi  =  Po  +  Pi*ii  +  -  +  Mik  i  =  l,2,...,n  (3-1) 


where 

n  is  the  number  of  sample  observations  in  our  data  set 

yi  is  the  i^^  observation  of  the  response 

Ej  is  the  model  error  associated  with 

Pq,  Pj,  ...,  p,^  are  the  unknown  linear  parameters. 


This  general  linear  model  can  be  expressed  in  matrix  form 

Y  =  XB  -l-E 


(3-2) 
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In  matrix  form,  the  least  squares  estimate  for  the  regression  coefficient  vector,  B,  is  giv¬ 
en  by  the  following  where  the  superscript  T”  indicates  matrix  transpose 

B  =  X’V.  (3-3) 

The  estimated  regression  equation  is  then 

Y  =  XB  (3-4) 


and  E  =  Y  -  Y  is  the  vector  of  model  residuals. 

Tb  determine  the  probabilities  of  cloud  amount  categories,  we  used  the  results  of 
the  general  linear  regression  model  described  above.  Our  goal  was  to  be  able  to  answer  the 
following  types  of  questions  that  are  critical  for  many  Air  Force  and  Army  mission  opera¬ 
tions: 

1.  Given  a  cloud  amount  forecast,  what  is  the  probability  that  the  observed  cloud 
amount  will  be  less  than  a  critical  threshold? 

2.  Given  a  cloud  amount  forecast,  what  is  the  probability  that  the  observed  cloud 
amount  will  be  within  a  critical  range  of  values? 

Td  answer  these  questions  we  must  understand  the  characteristics  of  the  error  dis¬ 
tribution  in  our  model,  lb  find  the  probability  that  the  observed  cloud  amount  will  be  less 

A  A  A 

than  or  equal  to  a  critical  threshold,  c  (i.e.,  Pr  |Y  ^  c|),  first  calculate  Y  =  XB  where  B 
is  the  best  estimate  for  the  regression  coefficients  using  least  squares  estimation.  Then, 
the  probability  can  be  computed  from  the  residual  distribution  as  follows: 

Pr  {  Y  :s  c  I  Y1 
=  Pr  { Y  +  E  ^  c  I  Y) 

=  Pr  {  E  <  c-Y  I  Y}. 

Figure  26  shows  graphically  how  to  compute  the  probability  in  Eq.  3-5  from  a  sample  error 
distribution;  the  shaded  area  under  the  curve  is  the  probability  referred  to  in  Eq.  3-5 
above.  It  is  a  trivial  extension  of  that  equation  to  compute  the  probability  that  the  ob¬ 
served  cloud  cover  is  vrithin  a  specified  range  by  finding  the  difference  between  the  proba¬ 
bilities  of  the  two  extreme  thresholds  that  define  the  range. 
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Linear  regression  assumes  that  model  residuals,  or  errors,  are  normally  distrib¬ 
uted  over  the  entire  sample:  E  ~  N(0,o^).  Our  results  have  verified  that  when  taken  as 
a  whole,  they  are  indeed  approximately  normal  (see  Figure  42  in  Section  3.4.1).  A  naive 
approach  to  computing  the  probabilistic  forecast  would  assume  that  the  residuals  are  nor¬ 
mally  distributed  across  any  range  of  cloud  cover  estimates.  However,  over  any  range  of 
estimates,  the  associated  residuals  may  be  far  from  normal  due  to  the  bounded  nature  of 
the  cloud  cover  variable.  In  most  cases  the  error  distributions  (especially  for  clear  and 
overcast  forecasts)  are  in  fact,  U-,  J-,  or  L-shaped.  (See  Figure  43  in  Section  3.4.2) 

Therefore,  a  better  way  to  compute  the  probability  in  Eq.  3-5  is  to  use  the  actual 
error  distributions  for  the  development  data  set.  We  divided  the  cloud  amount  forecasts 
into  11  bins:  0-4%,  5-14%, . .  . ,  85-94%,  95-100%  and  then  generated  histograms  of  the 
model  errors  for  each  of  the  forecast  categories.  The  probability  that  the  observed  cloud 
cover  was  less  than  a  critical  threshold,  given  the  estimated  cloud  coverage,  was  computed 
by  summing  the  error  histograms  corresponding  to  that  estimated  coverage. 

3.2  SPATIAL  CLUSTERING  OF  CLOUD  COVER 

The  first  step  in  the  model  development  sequence  was  to  reduce  the  number  of  re¬ 
gression  models  using  principal  components  analysis.  Regression  models  were  built  for 
each  cluster  rather  thai  Tor  individual  grid  points.  This  reduced  the  overall  number  of 
models  to  be  developed  and  also  increased  the  size  of  the  development  data  sets.  Cluster- 
based  development  did  not  reduce  the  specificity  of  the  final  models  since  parameters 
used  to  define  the  clusters  were  also  included  in  the  pool  of  potential  predictors. 

Clustering  was  constrained  to  individual  RTNEPH  boxes  and  used  location,  ter¬ 
rain,  and  cloud  amount  frequency  distribution  variables  for  each  of  the  4096  (i.e.,  64  x  64) 
RTNEPH  grid  points  in  a  box.  Specifically,  the  cluster  data  sets  included,  for  each  eighth- 
mesh  grid  point;  latitude,  longitude,  terrain  elevation,  terrain  type,  mean  cloud  cover,  and 
scale  distance.  The  last  two  variables  in  this  list  are  computed  from  the  RTNEPH  analy¬ 
ses  over  one  season  (3  months)  for  8  times  of  day. 

Scale  distance  is  a  parameter  of  the  Burger  Area  Algorithm  (BAA)  which  was  devel¬ 
oped  at  the  Air  Force  CJeophysics  Laboratory  and  is  described  in  Ref  5.  The  BAA  is  a  model 
of  the  probability  that  the  areal  cloud  coverage  is  less  than  or  equal  to  a  threshold  cover¬ 
age  value  given  the  area,  scale  distance,  and  mean  clearness.  The  scale  distance  can  be 
obtained  by  using  the  inverse  BAA  when  provided  with  a  cloud  cover  data  set.  It  is  a  mea¬ 
sure  of  the  range  of  significant  spatial  autocorrelation  of  cloud  cover. 
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We  performed  cluster  analysis  using  principal  components  for  four  box/season  com¬ 
binations  including  one  polar,  one  tropical,  and  two  midlatitude  boxes.  Principal  component 
analysis  (PCA)  has  two  main  applications;  1)  to  reduce  the  number  of  variables  in  a  data  set 
and  2)  to  detect  relationships  between  variables.  We  used  both  applications  in  our  analysis. 

Using  the  STATISTICAL*  factor  analysis  module,  we  extracted  principal  compo¬ 
nents  from  our  four  development  data  sets,  excluding  the  terrain  type  variable.  Succes¬ 
sive  principal  components  are  extracted  from  the  data.  The  first  component  is  a  best  fit 
through  the  data  that  maximizes  the  variance  of  that  component,  the  second  component 
is  orthogonal  to  the  first  and  maximizes  the  remaining  variability,  and  so  forth. 

Consecutive  principal  components  account  for  less  and  less  variability,  thus  it  is 
possible  to  reduce  the  number  of  variables  by  deciding  to  keep  only  the  first  few  compo¬ 
nents.  The  decision  to  cut  off  some  components  can  be  somewhat  subjective,  but  a  rule  of 
thumb  is  to  keep  those  components  which  account  for  at  least  as  much  variance  as  what 
one  individual  variable  would  account  for.  The  variance  of  the  individual  variables  (once 
normalized)  is  one.  Therefore  if  the  variance  of  a  principal  component  of  the  normalized 
data  set  is  greater  than  or  equal  to  one,  we  will  retain  it. 

Figures  27  and  28  show  the  variance  of  the  five  principal  components  calculated  for 
two  winter  cases:  box  38  in  the  midlatitudes  and  box  61  in  the  tropics  (refer  to  Figure  3 
for  the  RTNEPH  box  definitions).  The  first  figure,  corresponding  to  box  38,  shows  a  case 
in  which  the  first  two  principal  components  have  variance  greater  than  one.  We  retained 
both  components  for  further  analysis.  The  second  figure,  corresponding  to  box  61,  is  an 
example  of  the  subjectivity  involved  in  deciding  how  many  principal  components  to  retain 
for  analysis.  In  that  case  we  chose  to  retain  the  first  three  components  even  though  the 
variance  of  the  third  component  was  slightly  less  than  one. 

I*rincipal  components  are  linear  combinations  of  the  input  variables.  Figures  29 
and  30  show  the  normalized  breakdown  of  the  components  into  their  parts.  In  these  two 
cases,  the  first  principal  component  was  strongly  dependent  on  the  mean  cloud  cover  and 
latitude.  The  scale  distance  parameter  weighed  heavily  in  the  first  component  for  box  38 
and  in  the  second  component  for  box  61. 

Although  all  that  we  have  discussed  so  far  in  this  section  is  the  use  of  PCA  for  data 
reduction,  the  goal  of  this  clustering  task  was  not  to  reduce  the  five  terrain  and  cloud  dis¬ 
tribution  parameters  to  two  or  three  linear  combinations,  but  to  determine  clusters  of  grid 
points  with  similar  characteristics.  We  used  the  principal  components  to  determine  how 
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similar  grid  points  were  by  looking  for  groups  of  grid  points  in  a  scatterplot  of  the  principal 
components.  If  the  points  were  close  in  principal  component  space  they  were  also  close  in 
terms  of  cloud  climatology. 

We  found  that  the  grid  points  grouped  strongly  by  terrain  type  (land  and  water),  which 
makes  intuitive  sense  since  cloud  climatology  is  quite  distinct  in  water  and  land  regions  and 
the  bulk  of  the  grid  points  with  the  same  terrain  type  are  spatially  close.  Coastal  grid  points 
seemed  to  be  mixed  between  the  two  groups  as  mq)ected.  Figures  31  and  32  show  the  scatter- 
plots  for  boxes  38  and  61.  The  figures  for  box  38  are  two-dimensional  scatterplots  because 
we  retained  only  the  first  two  principal  components  in  our  analysis,  on  the  other  hand,  the 
scatterplots  for  box  61  are  three-dimensional  to  accommodate  the  three  components  we 
retained  for  analysis. 

The  two  other  cases  we  studied  showed  similar  results:  both  boxes  45  and  28 
showed  strong  clustering  by  terrain  type.  Although  not  conclusive,  this  analysis  led  us  to 
believe  that  we  could  generate  MOS  equations  for  distinct  terrain  types  within  RTNEPH 
boxes.  Ideally,  with  much  greater  resources,  it  would  be  best  to  develop  models  for  each 
grid  point  individually  and  only  after  comparing  the  models,  group  points  with  similar 
models  into  clusters.  However,  that  type  of  development  would  require  a  much  larger  de¬ 
velopment  data  set  and  modeling  budget. 

Findings  in  this  task  led  us  to  develop  a  number  of  terrain-specific  regression  mod¬ 
els  for  analysis.  However,  for  the  bulk  of  the  models  that  we  developed,  we  combined  land 
and  water  grid  points  to  provide  a  greater  number  of  data  points  for  model  development. 


3^  BUILDING  THE  POOL  OF  POTENTIAL  PREDICTORS 

The  selection  of  variables  to  use  as  potential  predictors  was  a  critical  task  in  the 
model  development  process.  We  selected  parameters  based  on  input  from  personnel  from 
Phillips  Laboratory,  the  Techniques  Development  Laboratory  (TDL)  of  the  National 
Weather  Service,  and  our  own  cloud  modeling  experts. 

Potential  predictors  were  taken  primarily  from  the  GSM  half-mesh  output  fields 
and  included  explanatory  variables  such  as  relative  humidity,  temperature,  and  surface 
pressure  as  well  as  derived  quantities,  transformed  variables,  and  differeni^  variables 
(difference  over  3-,  6-,  and  9-hour  periods).  Other  predictors  included  eighth-mesh  cloud 
amount  analyses  from  the  RTNEPH,  eighth-mesh  elevation  and  type  from  the  RTNEPH 
terrain  database,  latitude,  longitude,  date,  and  solar  zenith  angle. 
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Our  challenge  was  to  limit  the  total  number  of  variables  to  include  in  the  predictor 
pool  without  omitting  good  predictors  of  cloud  cover  at  eighth-mesh  resolution.  Limiting 
the  overall  number  of  variables  was  important  since  our  findings  may  lead  to  future  devel¬ 
opment  of  an  operational  cloud-amount  forecast  model  in  which  accuracy  and  economy  are 
both  major  concerns. 

Potential  predictors  included  point  variables  (e.g.,  surface  pressure  and  cloud 
amount),  column  variables  (e.g.,  temperature  and  wind  at  vertical  sigma  levels),  field 
variables  (e.g.,  terrain  elevation),  and  variables  at  times  other  than  the  valid  time.  Some 
of  the  data  fields  were  highly  correlated.  It  was  our  goal  to  reduce  the  number  of  potential 
predictors  by  taking  advantage  of  this  data  redundancy  using  layer  averaging,  principal 
components,  etc. 

The  set  of  potential  predictors  we  settled  on  included  approximately  eighty  vari¬ 
ables.  The  fifteen  most  significant  of  these  potential  predictors  were  selected  for  each  box/ 
season/time  using  a  “forward  selection”  scheme.  These  fifteen  were  then  used  as 
predictors  in  final  model  development.  Forward  selection  methods  (discussed  in  Ref  10) 
represent  one  way  to  reduce  the  overall  predictor  pool.  In  those  methods,  variables  are 
selected  one  at  a  time.  Once  a  variable  is  selected,  it  is  retained  without  chance  of  later 
being  discarded.  First  the  best  single  variable  is  selected.  With  k  potential  predictors,  this 
requires  k  simple  (one-variable)  regressions.  Next,  the  second-best  variable,  given  the 
particular  choice  of  the  best,  is  selected  by  performing  k-1  two-variable  regressions.  This 
process  is  continued,  adding  one  variable  at  each  step,  until  the  selection  of  an  additional 
variable  is  not  warranted  according  to  a  pre-specified  stopping  rule  or  until  all  variables 
have  been  selected. 

Our  stopping  rule  was  based  on  the  reduction  in  variance.  We  found  that  the  first 
five  to  ten  variables  accounted  for  the  majority  of  the  total  reduction  in  variance  achieved 
with  all  the  predictors.  After  looking  at  the  results  from  many  regions,  seasons,  and  times, 
we  decided  to  select  the  strongest  fifteen  variables  for  use  in  model  development.  By  keep¬ 
ing  fifteen  predictors,  we  accounted  for  the  bulk  of  the  reduction  in  variance  and  drastical¬ 
ly  reduced  the  computational  time  of  performing  the  regression.  In  addition,  we  retained 
enough  predictors  as  to  not  lose  critical  insights  into  the  subtle  differences  in  cloud  pro¬ 
cesses  in  varying  locations,  at  varying  times  of  day,  etc. 

In  Section  3.4,  we  present  our  modeling  results.  For  each  case,  we  show  which  fif¬ 
teen  predictors  were  selected  using  the  forward  selection  method  and  what  the  relative 
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ranking  of  those  variables  was.  It  would  probably  not  be  necessaiy  to  use  all  fifteen  predic¬ 
tors  in  an  operational  model,  but  the  increased  computations  for  this  feasibility  study 
were  well  worth  the  effort. 

A  complete  list  of  the  variables  used  in  this  study  is  included  in  Appendix  A.  We 
briefly  describe  the  RTNEPH  variables  below  and  then  devote  the  remainder  of  this  sec¬ 
tion  to  describing  the  GSM-based  predictors. 

3.3.1  Candidate  Predictors  From  the  RTNEPH  Cloud  Analyses  and 
Terrain  Database 

These  predictors  included  total  cloud  amounts  for  the  initial  time  and  for  times  at 
three-hour  intervals  preceding  the  initial  time  back  to  24  hours  prior  to  the  forecast  valid 
time.  The  MOS  equations  used  these  variables  to  model  persistence  effects  and  trends  in 
cloud  amount.  Terrain  elevation  and  type  were  also  included  in  the  pool  of  predictors  along 
with  the  terrain  gradient  in  the  direction  of  GSM  winds. 

A  key  predictor  based  on  the  RTNEPH  gridded  analyses  turned  out  to  be  advected 
cloud  cover.  Intuition  tells  us  that  to  forecast  cloud  cover  at  a  point  we  should  look  upwind 
to  see  what  weather  (i.e.,  cloud  cover)  will  move  into  the  area  around  that  point  over  the 
forecast  time  period. 

We  defined  a  relatively  simple  advected  cloud  cover  predictor  driven  by  forecast 
winds  at  the  valid  time  for  the  point  of  interest.  We  computed  average  low,  middle,  and 
high  layer  winds  (see  below  for  how  we  defined  the  three  vertical  regions)  over  a  horizon¬ 
tal  region  including  the  point  of  interest  and  the  eight  surrounding  half-mesh  grid  points. 
We  then  assumed  the  winds  acted  constantly  over  the  forecast  period  (3, 6,  or  9  hours)  to 
find  the  upwind  cloud  cover.  We  determined  the  upwind  cloud  amount  over  a  region  the 
size  of  a  half-mesh  grid  point  by  averaging  the  nearest  16  eighth-mesh  cloud  amount  val¬ 
ues.  We  produced  advected  cloud  amounts  for  low,  middle,  and  high  layer  winds  over  3-, 
6-,  and  9-hour  periods.  Further  discussion  of  these  predictors  and  their  frequency  of  occur¬ 
rence  in  the  strongest  fifteen  predictors  can  be  found  in  Section  3.4. 

3.3J2  Candidate  Predictors  From  the  Global  Spectral  Model 
Raw  Variables 

The  GSM  outputs  available  for  use  in  this  modeling  effort  included:  u  and  v  wind 
components  and  temperature  at  the  ten  lowest  vertical  (sigma)  levels,  relative  humidity 
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at  the  six  lowest  sigma  levels,  surface  pressure  and  model  surface  height  for  each  grid 
point  on  the  half-mesh  grid.  Output  variables  from  the  OZ-  and  12Z-initialized  model  runs 
were  both  available,  in  general,  though  frequently  data  from  only  one  model  run  were 
used  to  develop  the  MOS  equations  to  simulate  operational  constraints  (see  Section  2.2.2). 

Two  other  variables  were  obtained  directly  from  the  GSM  layer  data.  They  are;  the 
maximum  relative  humidity,  and  the  humidity  at  the  layer  above  the  maximum  (this  last 
variable  was  suggested  in  Ref.  20). 

Derived  Variables 


It  is  possible  to  derive  a  number  of  meteorologically  important  variables  from  the 
“raw”  u,  V,  temperature,  and  humidity  variables  mentioned  above.  An  obvious  predictor 
for  cloud  amount  that  can  be  derived  from  these  quantities  is  vertical  velocity.  Others  in¬ 
clude:  divergence,  vorticity,  temperature  advection,  layer  thickness,  wind  speed,  vertical 
wind  shear,  and  a  number  of  stability  parameters.  Each  of  these  derived  quantities  is  de¬ 
scribed  below. 


Divergence  —  is  a  measure  of  how  a  fluid  is  expanding  or  contracting.  In  this  ef¬ 
fort,  it  is  a  measure  of  iresoscale  motion  in  the  atmosphere.  In  a  simple  way,  convergence 
(negative  divergence)  can  be  associated  with  colliding  air  masses,  and  positive  vertical 
motion  which  are  both  indicative  of  increased  cloud  cover.  Positive  divergence,  on  the  oth¬ 
er  hand,  can  be  a  sign  of  clearing.  Horizontal  divergence  is  computed  as  the  sum  of  the 
partial  derivatives  of  the  u  and  v  wind  components,  where  u  corresponds  to  the  component 
of  the  wind  in  the  x-direction  (zonal  flow),  and  v  corresponds  to  the  component  of  the  wind 
in  the  y-direction  (meridional  flow) 


Divergence  =  • 


_  du  .  dv 
dx  dy 


(3-6) 


where  the  subscript  “h”  specifies  the  horizontal  components  only.  In  practice,  we  can  only 
approximate  this  theoretical  definition  by  computing  the  divergence  at  each  half-mesh 
grid  point  from  winds  surrounding  grid  points.  Figure  33  shows  the  geom#’''  ry  we  used 
to  implement  the  horizontal  derivatives.  The  divergence  at  the  center-point  'rf  Jie  grid  is 
given  by 


DiVh 


(Ua  -  Uc)  .  ~ 

2Ax  2Ay 


(3-7) 


where  Ax  and  Ay  are  the  distances  between  adjacent  grid  points  in  the  x  and  y  directions, 
respectively.  Divergence  is  very  sensitive  to  errors  in  calculation  because  the  two  terms 


51 


in  Eq.  3-7  are  usually  small,  often  of  opposite  sign,  and  often  very  close  in  magnitude. 
Thus  we  need  to  know  the  terms  very  accurately  in  order  to  make  a  reasonable  estimate 
of  the  divergence. 

Vorticity  —  is  a  measure  of  the  local  rotation  of  a  fluid.  It  is  intimately  connected 
to  divergence  because  as  fluids  converge  or  diverge  their  rate  of  rotation  changes  to  con¬ 
serve  angular  momentum.  Mathematically,  vorticity  is  given  by 

Voracity  =  k  .  (V,  X  7,)  =  g  -  (3.8) 


Implementation  used  approximate  derivatives  to  determine  the  vorticity  at  the  center  of 
the  grid  (refer  again  to  Figure  33): 


Vorth  = 


(Va  -  Vc) 
2Ax 


<Ud  ~  Ub> 
2Ay 


(3-9) 


Vertical  Velocity  —  is  a  measure  of  the  vertical  motion  in  the  atmosphere.  It  is 
usually  a  good  measure  of  convective  activity  and  precipitation  and  thus,  indirectly,  doud 
cover.  It  is  usually  derived  from  divergence,  and  thus  is  sensitive  to  errors  in  the  diver¬ 
gence  computation.  There  are  many  methods  to  compute  vertical  velocity.  We  used  the 
method  followed  in  the  Air  Force  version  of  the  Global  Spectral  Model  (Ref  2).  The  vertical 

velocity,  6 ,  is  the  total  derivative  of  the  sigma  coordinate  as  follows; 

d  =  o(C  -t-  D)  -  (C^  If)  (3-10) 

where 


o 


0 


C  =  V  •  Vq,  q  =  log(pressure) 
D  =  Divergence 

and  where  d  =  0  at  o  =  0  and  a  =  1. 
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In  our  implementation  of  this  equation,  we  limited  the  vertical  range  over  which 
we  performed  the  integrations  to  the  ten  reported  sigma  levels  available  to  us.  That  is,  the 
part  of  the  atmosphere  above  the  top  sigma  level  of  our  data  (approximately  120  mb)  was 
ignored.  Due  to  the  nature  of  the  MOS  approach,  this  approximation  was  most  likely  of 
little  consequence,  because  MOS  can  use  the  trend  in  vertical  velocity  rather  than  the  ab¬ 
solute  value  as  a  predictor  of  cloud  amount. 

Temperature  Advection  —  is  a  measure  of  how  fast  the  temperature  gradient  is 
changing.  It  can  be  a  measure  of  baroclinic  development  and  thus  of  the  changes  in  cloud 
cover.  Temperature  advection  is  given  by  the  following  equation: 

Temperature  Advection  •  Vj^T  =  (3-11) 

Again,  we  approximate  the  partial  derivatives  using  the  gridded  data  available  to  us  and 
implement  temperature  advection  at  the  center  of  the  grid  using  the  following  (refer  again 
to  Figure  33): 

laaVh  -  Uavg  2Ax  2Ay  (3-12) 

where  Uavg  and  Vavg  are  the  u  and  v  wind  components  averaged  over  the  center  point  o  and 
its  two  neighbors  along  each  direction  (we  weight  the  center  point  twice  as  heavily  as  the 
other  grid  points). 

Wind  Speed  —  is  simply  the  resultant  speed  computed  using  the  u  and  v  wind 
components  at  the  surface  as  follows: 

Wind  Speed  =  (3-13) 

Vertical  Wind  Shear  —  is  a  measure  of  the  magnitude  of  the  wind  shear  between 
high  and  low  layers  in  the  atmosphere.  We  computed  wind  shear  as 

Wind  Shear  =  y(Uhigh  -  +  (v^gh  "  (3-14) 


where  the  high  layer  winds  are  averages  over  the  top  five  sigma  levels  (approximately  375 
to  100  mb)  and  the  low  layer  winds  are  averages  of  the  lower  two  sigma  levels  (surface  to 
approximately  800  mb). 
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Stability  Measure  —  is  a  simple  measure  of  the  temperature  change  with  height 
in  the  atmosphere.  We  computed  two  stability  measures  as  follows: 

StabiUty  Measure  1  =  (3.15) 

Stability  Measure  2  =  Tgsomb  (3-16) 


where  all  pressures  are  approximate  because  the  temperatures  are  actually  evaluated  on 
sigma  coordinates  which  are  defined  relative  to  surface  pressure. 

Layer  Thickness  —  is  the  geometric  thickness  of  the  layer  bounded  at  the  bottom 
by  pressure  Pb  and  at  the  top  by  pressure  Pt.  It  is  a  type  of  stability  measure  in  the  atmos¬ 
phere.  The  algorithm  to  compute  the  thickness  involves  numerical  integration  of  the  hy¬ 
drostatic  equation.  It  is  described  in  Ref.  25.  We  used  Pb=P8urface  and  Pt=850  mb. 

IVansformed  Variables 

Another  group  of  potential  predictors  included  transformed  variables.  During  pre¬ 
liminary  data  analysis  we  studied  scatterplots  of  individual  predictors  versus  cloud 
amount  and  residual  plots  which  led  us  to  try  non-linear  transformations  of  key  variables, 
especially  relative  humidity  (RH).  Through  trial  and  error  and  analysis  of  previous  re¬ 
sults  we  settled  on  the  following  group  of  transformed  predictors: 

•  log(RH) 

•  RH2 

•  RH^ 

•  vertical  velocity  x  RH 

•  RH  >  50%,  RH  >  70%,  RH  >  90%. 

The  f  rst  three  of  these  predictors  emphasize  extreme  relative  humidity  values.  In 
the  same  way,  the  product  of  vertical  velocity  and  relative  humidity  emphasizes  those  sit¬ 
uations  when  both  variables  are  extremely  low  or  high.  This  predictor  was  suggested  by 
personnel  at  the  TDL  (Ref.  11).  The  same  group  has  also  found  thresholded  predictors  use¬ 
ful  in  building  MOS  equations  using  output  fh>m  the  Nested  Grid  Model.  The  last  predic¬ 
tors  are  examples  of  thresholded  or  binaiy  variables,  equal  to  one  if  the  humidity  is 
greater  than  the  threshold,  and  zero  otherwise. 
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Differenced  Variables 

The  Global  Spectral  Model  generates  forecasts  every  twelve  hours  in  three  hour  in¬ 
crements.  In  some  cases  it  may  be  desirable  to  use  model  output  from  three  hours  before 
or  after  the  actual  valid  time  in  the  predictor  pool  to  account  for  possible  spin-up  problems 
in  the  model.  Due  to  computational  limitations  we  were  not  able  to  add  data  from  both 
before  and  after  the  valid  time  in  our  pool  of  potential  predictors.  Instead,  we  used  3-,  6-, 
and  9-hour  differences  in  key  variables  to  capture  the  temporal  component  of  the  model. 
Using  time  differenced  predictors  in  the  MOS  equations  allowed  us  to  avoid  any  absolute 
errors  in  the  model  forecasts  (e.g.,  moisture  values  are  known  to  be  biased  at  high  sigma 
levels,  Ref  20)  and  instead  rely  on  trends  in  meteorological  variables.  For  example,  when 
making  a  6-hour  forecast  valid  at  15Z,  we  use  the  GSM  forecast  variables  valid  at  15Z  and 
the  change  in  forecast  variables  over  the  9Z  to  15Z  period.  The  difference  predictors  can 
be  used  to  predict  changes  in  cloud  amount  from  the  amount  at  initial  time.  This  is  the 
strength  of  our  MOS  approach:  combining  observations  at  initial  time,  forecast  variables 
at  the  valid  time,  and  differences  in  model  variables  from  initial  to  valid  time,  so  as  to  ac¬ 
commodate  potential  model  biases  and  spin-up  problems. 

Reduction  of  the  Number  of  Variables 

If  we  considered  all  of  the  raw,  derived,  transformed,  and  differenced  variables 
mentioned  above  for  each  vertical  level,  our  pool  of  potential  predictors  would  contain  over 
200  variables  for  each  grid  point.  Because  computations  involved  in  forward  variable 
selection  and  regression  in  general  are  so  costly,  we  investigated  the  tradeoffs  associated 
with  reducing  the  overall  number  of  potential  predictors. 

Based  on  our  preliminary  analyses,  we  concluded  that  we  could  reduce  the  ten  sig¬ 
ma  layers  to  averaged  low,  middle,  and  high  layers  with  little  loss  of  fidelity.  In  fact,  we 
found  that  the  three-layer  variables  were  more  stable  predictors  of  total  cloud  amount 
than  their  ten  individual  layer  counterparts.  This  enabled  us  to  cut  by  over  two-thirds  the 
number  of  column  variables  in  the  pool  of  potential  predictors.  The  definition  of  the  low, 
middle,  and  high  layers  we  used  in  this  effort  with  respect  to  the  ten  GSM  sigma  layers 
is  shown  in  Figure  34. 

The  resulting  list  of  potential  predictors  used  in  our  final  regression  model  develop¬ 
ment  is  found  in  Appendix  A. 
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3.4  REGRESSION  MODEL  RESULTS 


Modeling  consisted  of  two  phases:  exploratory  and  final.  In  the  exploratory  phase 
discussed  in  previous  sections,  we  built  up  the  predictor  list,  studied  the  effects  of  time  of 
day,  length  of  forecast,  and  location  on  data  availability,  and  defined  the  general  modeling 
process.  During  that  phase  we  mainly  studied  three  regions:  RTNEPH  boxes  44,  45,  and 
61.  Most  of  our  results  from  that  phase  were  qualitative. 

In  the  final  modeling  phase  we  chose  a  number  of  box/season/time  combinations 
based  on  availability  of  timely  data,  variety  of  terrain  types  and  cloud  climatologies,  and 
strategic  importance.  A  base  set  of  59  combinations  was  selected  for  analysis.  Td  that  set 
we  added  12  more  cases  that  were  terrain  specific;  that  is,  they  only  included  grid  points 
of  a  single  terrain  type:  land  or  water  as  suggested  by  our  cluster  analysis.  All  71  analysis 
cases  are  listed  in  Table  3. 

We  defined  a  standard  set  of  predictors  (listed  in  Appendix  A)  and  i  odeling  opera¬ 
tions  including: 

•  selecting  the  strongest  predictors  from  the  pool 

•  determining  regression  coefficients 

•  computing  partial  correlations 

•  determining  the  reduction  in  variance 

•  analyzing  model  residuals. 

From  model  results  we  then  determined  categorical  cloud  amount  probabilities.  Finally, 
we  compared  the  MOS  results  to  persistence  using  meteorological  skill  scores  and  pre¬ 
dicted  future  model  performance. 

Regression  models  were  developed  from  data  sets  containing  one  season  of  data 
and  225  eighth-mesh  grid  points  (fewer  for  terrain-specific  models)  one  grid  point  for  each 
of  15  X  15  half-mesh  cells.  (The  original  GSM  outputs  were  on  17  x  17  half-mesh  grids  af¬ 
ter  processing  at  Phillips  Laboratory.  We  used  only  the  inner  15  x  15  cells  so  as  to  be  able 
to  compute  the  horizontal  gradients  used  in  divergence,  vorticity,  and  temperature  advec- 
tion  computations).  The  four  seasons  were  defined  as  follows: 

Winter  =  (December,  January,  February} 

Spring  =  (March,  April,  May} 

Summer  =  (June,  July,  August} 

Fall  =  (September,  October,  November}. 
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Table  3  List  of  7 1  Model  Development  Data  Sets 


CASE* 

1 


8 

9 

"io 

11 

IF 

"U 

14 

15 

"is 

"iT 

18 

"19 

20 

21 

22 

24 

25 

26 

Is 

"m 

30 

33 

34 

35 

TF 

39 

40 


CASE  NAME 


30/SP/3/3/12 


30/SP/3/6/12 


30/SP/3/9/12 


30/SP/9/3/0 


BOX  SEASON 


SPRING 


TIME  OF 
DAY  (2) 


FORECAST 

LENGTH 

(HOURS) 


GSM 

INITIAL 


TERRAIN 

TYPE 


SPRING 


SPRING 


SPRING 


30/SP/9/6/0 

30 

SPRING 

30/SP/9/6/12 

30 

SPRING 

30/SP/9/9/12 

30 

SPRING 

30/SP/21/3/12 

30 

SPRING 

30/SP/21/6/12 

30 

SPRING 

30/SP/21/6/0 

30 

SPRING 

30/SP/21/9/0 

30 

SPRING 

21 

9 

0 

ALL 

30/SU/9/3/0 

30 

SUMMER 

9 

3 

0 

ALL 

30/SU/9«/0 

30 

SUMMER 

9 

6 

0 

ALL 

30/SU/9/6/12 

30 

SUMMER 

9 

6 

12 

ALL 

30/SU/9/9/12 

30 

SUMMER 

9 

9 

12 

ALL 

30/SU/21/3/12 

30 

SUMMER 

21 

30/SU/21/6/12 

30 

SUMMER 

21 

30/SU/21/6/0 

30 

SUMMER 

21 

30/SU/21/9/0 

30 

SUMMER 

21 

30/FA«/3/0 

30 

FALL 

9 

30/FA/9/6/0 

30 

FALL 

9 

30/FA/9/6/12 

30 

FALL 

9 

30/FA/9/9/12 

30 

FALL 

9 

30/FA/21/3/12 

30 

FALL 

21 

30/FA/21/6/12 

30 

FALL 

21 

30/FA/2 1/6/0 

30 

FALL 

21 

30/FA/2 1/9/0 

30 

FALL 

21 

30/WI/9/3/0 

30 

WINTER 

9 

30AA/I/9/6/0 

30 

WINTER 

9 

30/WI/9/6/12 

30 

WINTER 

9 

30AWI/9/9/12 

30/WI/21/3/12 

30/WI/21/6/12 

30/WI/21/6/0 

30AWI/21/9/0 

44/SU/0/3/12 

44/SU/0/6/12 

44/SU/0/9/12 

44/SU/1 5/3/0 

44/SU/1 5/6/0 
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CASE  NAME 

BOX 

SEASON 

44/SU/1 5/9/0 

44 

SUMMER 

44/WI/0/3/12 

44 

WINTER 

44/WI/0/6/12 

44 

WINTER 

44/WI/0/9/12 

44 

WINTER 

44/WI/1 5/3/0 

44 

WINTER 

44/WI/1 5/6/0 

44 

WINTER 

44/WI/1 5/9/0 

44 

WINTER 

61/SU/3«/12 

61 

SUMMER 

61/SU/3/6/12 

61 

SUMMER 

61/SU/3/9/12 

61 

SUMMER 

61/SU/1 5/3/0 

61 

SUMMER 

61/SU/1 5/6/0 

61 

SUMMER 

61/SU/1 5/9/0 

61 

SUMMER 

61/WI/3/3/12 

61 

WINTER 

61/WI/3/6/12 

61 

WINTER 

61/WI/3/9/12 

61 

WINTER 

61/WI/15/3/0 

61 

WINTER 

61/WI/15/6/0 

61 

WINTER 

61/WI/1 5/9/0 

61 

WINTER 

3Q/SU/21/3/12L 

30 

SUMMER 

30/SU/21/6/0L 

30 

SUMMER 

30/SU/21/9/0L 

30 

SUMMER 

30/SU/21/3/12W 

30 

SUMMER 

30/SU/21/6/0W 

30 

SUMMER 

30/SU/21/9/0W 

30 

SUMMER 

61/WI/3/3/12L 

61 

WINTER 

61/WI/3/6/12L 

61 

WINTER 

61/WI/3/9/12L 

61 

WINTER 

61/WI/3/3/12W 

61 

WINTER 

61/WI/3/6/12W 

61 

WINTER 

61/WI/3/9/12W 

61 

WINTER 

TMEOF 
DAY  (2) 


FORECAST 

LENGTH 

(HOURS) 


GSM 

INITIAL 

TtME(Z) 


TERRAIN 

TYPE 


ALL 


ALL 


ALL 


ALL 


ALL 


ALL 


ALL 


ALL 


ALL 


ALL 


ALL 


ALL 


ALL 


ALL 


ALL 


ALL 


ALL 


ALL 


ALL 


LAND 


LAND 


LAND 


WATER 


WATER 


WATER 


LAND 


LAND 


LAND 


WATER 


WATER 


WATER 


3.4.1  Predicting  Cloud  Cover  Fraction 


Strongest  Predictors 

We  generated  regression  models  for  all  71  cases  listed  in  Table  3.  The  first  step  in 
the  process  was  to  determine  the  strongest  predictors  in  each  case.  The  fifteen  strongest 
predictors,  as  selected  by  the  forward  stepwise  regression  module  in  the  STATISTIC  A™ 
software,  are  listed  for  each  of  the  71  model  cases  in  Table  4.  Brief  descriptions  of  the  vari¬ 
ables,  listed  by  their  abbreviated  names  in  the  table,  can  be  found  in  Appendix  A. 
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Table  4  List  of  Strongest  Fifteen  Predictors  Selected  Using 
Forward  Stepwise  Regression  for  the  71  Development  Data  Sets 
(Coded  Variable  Names  Listed  Here  are  Described  in  ^pendix  A) 
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1 

f 

> 

4 

B 

B 

» 

B 

44«U1S«0 

ACCMB 

RHMK 

caa 

8BCH 

ACCU 

NEPH24 

ACCHB 

DBVTH 

4«SUt&ttiO 

RHMX 

accbm 

ACCU 

8FCH 

NEPH24 

OBML 

TOIF1 

44MVMM3/12 

ACCL3 

NCPH3 

ACCH3 

ACCMB 

TU 

VTL 

44/WM)l«/12 

ACCU 

ACCMB 

NEPHB 

Qjgi 

TM 

RH70 

OVH 

NCPH24 

44MMM12 

ACCU 

cai 

caa 

ACCMB 

VTi 

TM 

NHTO 

OVH 

4«WV1S/M> 

ACCU 

csa 

ACCMB 

NEPHB 

TAN 

NCPHIB 

TM 

umiMW 

ACCU 

ggi 

ACCMB 

QB 

TM 

NCPHt2 

TAH 

44AW15/W0 

ACCU 

TM 

ACCMB 

TAH 

NEPHIB 

LON 

6t;suaa/l2 

ACCM3 

ACCH3 

LAT 

ACCU 

NEPH21 

NEPHIS 

NEPHB 

6i/suraMn2 

ACCMe 

UIRHM 

NEPN1B 

ACCU 

ACCHB 

NEPHIS 

TAM 

ZENTTH 

61/SU3^2 

ACCU 

sa 

ACCHB 

NCPHIB 

NEPH12 

ACCMB 

ca 

TYPt 

oi/surtsoo 

ACCM3 

lAT 

ACCM3 

ACCU 

MM 

TAM 

va 

61/SU1SW) 

ACCH6 

LAT 

NEPHa4 

ACTMO 

VMM 

VTL 

ESI 

omiH 

6i/sunsm> 

ACCM» 

NEPH24 

UT 

ACCHB 

EH 

VTL 

ESI 

OBMH 

61MW3fln2 

ACCMS 

NEPH21 

NEWS 

UT 

ACCU 

VTM 

ACCHB 

NHM2 

61MW3A/12 

ACCMB 

^22 

lAT 

ACCHB 

ACCU 

06VM 

VM 

NEPHB 

t^mvy9n2 

NEPH21 

ACCMB 

lAT 

ACCLB 

ACCHB 

OBVM 

VM 

DELEV 

«1AW1S/aO 

ACCMS 

NEPH24 

U 

ACCU 

ACCHB 

NEPHB 

NEPHB 

ACCMB 

MM 

06VM 

ACCHB 

ACCU 

NCPH12 

VTM 

61/WV1SW0 

ACCMB 

NEPH24 

ACCU 

VTM 

NEPH21 

RHL4 

OBOVL 

ACCHB 

30aUQ1/3n2L 

ACCM3 

|nCPH3 

ACCH3 

ELEV 

DVH 

NCPKBt 

DBP 

3»8U91«0(. 

ACCMB 

ESI 

NEBHB 

ELEV 

ACCHB 

OVH 

06VTH 

OVL 

3(VSUf21/ML 

ACCMB 

THICK 

NEPHB 

ACCHB 

VL 

ZCMTH 

VTL 

3QW21/V12W 

ACCM3 

ACCHB 

OVH 

03P 

OBVM 

GilH 

TAM 

3aSUI3l>MW 

ACCMB 

OVN 

yjjl 

ACCHB 

UM 

NCPHIB 

VTL 

DBVM 

3VSU3l/g^ 

ACCMB 

OVH 

AOCH9 

VTL 

NEMtB 

Q2I 

UM 

DAY 

61/vm3n2L 

ACCM3 

NEPH3 

NEPH21 

NHM 

ACCU 

ACCHB 

UT 

VTM 

61/WV3«n2L 

ACCMB 

NCPK21 

ACCHB 

UT 

ACCLB 

VTM 

LON 

NEPHB 

ei/WV3W12L 

ACCMB 

NEPH21 

ACCHB 

UT 

ACCLB 

LON 

DBWL 

OELEV 

61/WV3^2W 

ACCU 

NEPN21 

UL 

NEPHIS 

NEPHB 

ACCHB 

OBWH 

\NRHL 

81AW3«'12W 

NEPH2I 

ACCU 

UT 

NEPHIS 

WH 

VM 

ZENITH 

ACCHB 

61/WV3W12W 

NEPH15 

ACCU 

NEPH21 

UT 

WH 

VM 

VTM 

DBVL 
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Figure  35  shows  a  histogram  of  those  predictors  that  were  among  the  fifteen  stron¬ 
gest  most  often  over  all  71  development  data  sets.  In  that  figure  we  normalized  the  number 
of  occiurences  to  the  total  number  of  possible  occiirrences  for  each  variable.  RecaU  that  some 
predictors  (e.g.,  wind  speed)  are  common  to  all  3-,  6-,  and  9-hour  forecasts,  whereas  others 
(e.g.,  6-hour  change  in  surface  pressure)  are  unique  to  only  one  forecast  length. 

While  Figure  35  reflects  the  frequency  with  which  specific  predictors  appear  in  the 
strongest  fifteen,  it  does  not  contain  information  about  whether  the  variables  were  very 
strong  predictors  of  cloud  cover  (occurring  always  in  the  top  half  of  the  group)  or  simply  sup¬ 
porting  predictors  (occurring  in  the  lower  half  of  the  group  of  fifteen).  For  this  reason  it  was 
necessary  to  anal3rze  this  figure  in  conjunction  with  the  ordered  predictor  lists  of  Table  4  to 
make  more  specific  conclusions  about  the  meaning  and  utility  of  various  predictors. 

It  is  important  to  note  however,  that  especially  for  those  predictors  that  occur  in 
supporting  roles,  the  absolute  utility  of  any  individual  predictor  cannot  be  ascertained  by 
its  rank  position  in  the  top  predictors,  or  even  by  its  contribution  to  the  overall  reduction 
in  variance.  This  is  due  to  the  high  correlation  between  some  predictors.  The  forward  re¬ 
gression  technique  finds  the  set  of  fifteen  predictors  that  together  reduce  the  variance  by 
the  greatest  amount.  There  may  be  alternate  sets  of  predictors  that  perform  nearly  as  well 
with  different  supporting  predictors. 

Similarly,  for  a  different  development  data  set  we  may  find  that  the  group  of  sup¬ 
porting  predictors  changes.  We  performed  a  limited  analysis  in  which  we  split  a  develop¬ 
ment  data  set  in  half  and  found  the  strongest  predictors  from  the  two  halves  individually 
and  the  whole  data  set.  We  found  the  very  strongest  predictors  remained  constant  across 
all  three  data  sets.  Some  of  the  remaining  variables  changes  rank  position  or  dropped  out 
altogether  to  be  replaced  by  others. 

Keeping  these  caveats  in  mind,  and  realizing  that  our  conclusions  are  based  only 
on  one  year  of  data,  a  number  of  interesting  observations  drawn  from  Figure  35  and 
Table  4  along  with  results  from  our  exploratory  analysis  phase  are  grouped  loosely  by  pre¬ 
dictor  type  and  listed  below: 

Advected  Cloud  Cover  Terms 

•  The  advected  cloud  cover  terms  (ACCL3,  ACCM3,  ACCH3,  ACCL6,  etc.)  were 
very  strong  predictors  (frequently  in  the  strongest  five)  for  all  boxes,  times  of 
day,  seasons,  and  forecast  lengths  we  studied. 

•  The  advection  terms  had  their  greatest  influence  in  6-hour  forecasts,  especially 
in  box  30. 
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•  The  high-layer  advection  terms  did  not  perform  as  strongly  as  low-  and  middle- 

layer  advection  in  all  box  regions.  That  is,  although  they  were  frequently  in  the 
group  of  the  strongest  predictors,  they  often  played  supporting  roles. 

Persisted  Cloud  Cover  Amount 


•  The  cloud  cover  at  initial  time  (i.e.,  the  persistence  forecast)  was,  as  expected, 
a  strong  predictor  along  with  the  advection  terms  in  boxes  30  and  44. 

•  The  cloud  cover  at  initial  time  was  rarely  among  the  very  strongest  predictors 
for  box  61,  although  it  frequently  turned  up  in  a  supporting  role. 

•  The  cloud  cover  at  initial  time  generally  dropped  in  significance  with  forecast 
length  for  any  given  box/season/time  combination. 

•  The  24-hour  persisted  cloud  cover  turned  up  more  often  in  the  very  strongest 
predictors  for  box  61  than  in  boxes  30  and  44.  However,  24-hour  persistence 
occurred  more  often  overall  in  box  44. 

•  Cloud  cover  21  hours  before  valid  time  was  a  frequent  contributor  in  box  30, 
probably  representing  24-hour  persistence  that  was  three  hours  old  at  21 
hours  before  forecast  valid  time. 

•  In  box  61,  the  24-hour  persisted  cloud  cover  was  more  prevalent  at  15Z 
(approximately  10  a.m.  local  time)  than  at  3Z  (approximately  10  p.m.). 

Relative  Humidity  Variables 


•  Humidity  terms  in  general  were  very  significant  predictors. 

•  Maximum  layer  relative  humidity  occurred  frequently  in  the  strongest  fifteen 
for  all  box/season/time  combinations. 

•  Mid-level  relative  humidity  and  transformations  of  that  variable  (In(RH), 
RH^,  and  RH^)  were  also  present  frequently. 

•  Binary  relative  humidity  variables  were  also  in  the  strongest  predictors  for 
box  30. 

•  3-,  6-,  and  9-hour  humidity  change  variables  occurred  with  some  regularity  in 
box  30. 

Other  Variables 

•  Change  in  elevation  along  the  local  wind  direction  was  a  common  predictor. 

•  Latitude  was  an  important  predictor  for  box  61,  but  less  so  for  boxes  30  and  44. 

•  Mid-level  u  wind  and  low-level  v  wind  components  were  common  for  box  30. 
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•  3-,  6*,  and  9-hour  changes  in  surface  pressure  were  significant  predictors  in 
box  30. 

•  Low-level  vorticity  was  a  commonly-selected  predictor  across  all  three  regions. 

•  Mid-level  temperature  occurred  with  some  frequency  in  box  44. 

•  Mid-  and  high-level  temperature  advection  turned  up  frequently  in  box  44. 

•  Vertical  velocity  was  a  stronger  predictor  in  box  61  than  in  boxes  30  or  44. 

•  In  box  61,  the  vertical  velocity  was  more  prevalent  at  15Z  (approximately 
10  a.m.  local  time)  than  at  3Z  (approximately  10  p.m.). 

•  Elevation  and  stability  parameters  occurred  more  often  in  the  strongest  fifteen 
predictors  for  box  44  at  OZ  (approximately  6  p.m.)  than  at  15Z  (approximately 
9  a.m.),  characteristic  of  afternoon  heating  and  lifting  over  elevated  terrain. 

The  strength  and  ubiquitousness  of  the  simple  cloud  cover  advection  predictors 
surprised  us.  There  were  two  possibilities  for  their  predictive  strength:  1)  they  repre¬ 
sented  regional  average  cloud  cover  amounts  (these  variables  were  averaged  to  half-mesh 
resolution)  or  2)  they  really  represented  advection  of  regional  cloud  systems.  lb  determine 
whether  their  predictive  effects  were  due  to  advection  or  to  area  averaging  we  substituted 
winds  rotated  by  180°  in  the  computations.  The  result  was  that  the  advection  terms 
dropped  out  dramatically  from  the  strongest  predictors.  Not  only  did  this  imply  that  these 
predictors  represented  advection,  but  it  also  reflected  the  accuracy  of  the  GSM  winds. 

The  last  twelve  cases  in  Table  4  correspond  to  models  developed  for  specific  terrain 
types  (i.e.,  land  or  water).  The  strongest  predictors  for  these  terrain-specific  cases  were 
distinct  from  those  developed  for  the  same  regions  and  times  including  all  terrain  types. 
This  conclusion  supports  our  findings  during  the  clustering  task;  grid  points  of  similar 
terrain  type  which  grouped  tightly  in  principal  component  space  represented  unique  mod¬ 
el  development  cases  with  unique  cloud  distribution  characteristics.  (We  did  not  develop 
models  for  coastal  grid  points  because  there  were  not  enough  data  points  to  support  mean¬ 
ingful  results.) 

Many  of  the  same  observations  made  previously  and  discussed  above  for  the  59  all¬ 
terrain-type  cases  are  also  true  for  these  twelve  terrain-specific  cases  (e.g.,  the  overall 
strength  of  the  advection  and  persistence  terms).  However,  a  few  observations  specific  to 
these  cases  can  also  be  made: 

•  Box  61  land/water  cases  differed  less  than  box  30  land/water  cases  from  their 
all-terrain-t3q)e  counterparts. 

•  We  saw  a  greater  dependence  on  wind  and  divergence  predictors  in  the  water 
model  than  in  the  land  model  for  all  forecast  lengths  in  box  30/summer. 


•  Box  61/winter  showed  a  similar  emphasis  on  wind-derived  predictors. 

•  Cloud  cover  at  the  initial  time  was  a  stronger  predictor  over  land  than  over  wa¬ 
ter  in  box  30/summer. 

•  Temperature  and  temperature  advecdon  terms  occurred  in  all  three  water 
cases  in  box  30/summer. 

Regression  Coefficients 

After  determining  the  strongest  fifteen  predictors  for  each  data  set,  we  determined 
the  model  coefficients  (or  weights)  using  the  standard  linear  regression  technique  offered 
in  STATISTICA”* .  In  Tables  5  through  9  we  show  the  raw  regression  weights  (B)  and  the 
standardized  regression  weights  (BSTD)  for  5  of  the  71  models  (cases  24,  25,  27,  42,  and 
52  in  Table  3).  We  selected  these  five  cases  to  highlight  in  this  results  section  because  they 
represented  a  good  cross  section  of  the  data.  Three  of  the  cases  are  based  on  one  box/sea¬ 
son/time  combination,  but  have  different  forecast  lengths  (3, 6,  and  9  hours).  One  can  de¬ 
tect  the  effects  of  forecast  length  on  model  development  by  comparing  their  results.  The 
other  two  cases  are  from  different  boxes,  seasons,  and  times. 

The  magnitude  of  the  raw  B  weights  is  not  a  measure  of  the  relative  strength  of  the 
predictors,  but  rather  a  reflection  of  the  predictors*  units  and  magnitude.  The  standard¬ 
ized  regression  weights  are  those  we  would  have  obtained  had  we  first  standardized  all 
the  predictors  to  mean  0  and  standard  deviation  1.  Thus,  the  standardized  values  give  us 
a  measure  of  the  relative  contribution  of  each  predictor  to  the  overall  cloud  amount  predic¬ 
tion  equation.  Positive  weights  reflect  the  positive  partial  correlations  between  the  pre¬ 
dictor  and  the  predictand;  e.g.,  mid-level  relative  humidity  (RHM)  and  cloud  cover  in 
Table  7.  Negative  weights  imply  a  negative  relationship;  e.g.,  3-hour  change  in  surface 
pressure  (D3P)  and  cloud  cover  in  Table  5. 


64 


Table  5  Regression  Results  for  Case  #24  (Box  30,  Fall,  21Z,  3-Hour 
Forecast)  Showing  Top  15  Predictors,  Regression  Weights 
(Raw  and  Standardized),  and  Partial  Correlations 


VARIABLE  NAME 

B 

BSTD 

PARTIAL 

CORRELATION 

ACCL3 

.171 

.140 

7.62e-2 

NEPH3 

.307 

.300 

.252 

ACCH3 

.151 

.123 

7.93e-2 

NEPH24 

4.26e-2 

4.26e-2 

5.65e-2 

ZENITH 

-.161 

-5.29e-2 

-8.10e-2 

RHM 

3.48e-2 

1.71e-2 

1.61  e-2 

NEPH6 

6.08e-2 

5.61  e-2 

5.57e-2 

ACCM3 

.157 

.128 

5.82e-2 

D3RHM 

.331 

5.00e-2 

6.85e-2 

D3P 

-2.04e-4 

-4.01  e-2 

-5.77e-2 

RHL4 

.104 

5.02e-2 

-5.39e-2 

D3VL 

-1.19e-2 

-3.24e-2 

-5.08e-2 

DELEV 

4.36 

2.70e-2 

4.30e-2 

LON 

1.27e-3 

2.75e-2 

4.190-2 

NEPH15 

3.62e-2 

3.59e-2 

4.20e-2 
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Table  6  Regression  Results  for  Case  #25  (Box  30,  Fall,  21Z,  6-Hour 
Forecast)  Showing  Tbp  15  Predictors,  Regression  Weights 
(Raw  and  Standardized),  and  Partial  Correlations 


VARIABLE  NAME 

B 

BSTD 

PARTIAL 

CORRELATION 

ACCL6 

.248 

.193 

.127 

RHM 

8.74e-2 

4.296-2 

3.166-2 

NEPH6 

.210 

.194 

.172 

ACCH6 

.120 

9.256-2 

7.316-2 

ZENITH 

-.154 

-4.986-2 

-6.906-2 

NEPH21 

7.05e-2 

7.066-2 

8.606-2 

D6RHM 

.305 

7.686-2 

8.546-2 

RH70 

6.21  e-2 

7.056-2 

6.566-2 

D6P 

-1 .446-4 

-5.246-2 

-6.526-2 

ACCM6 

6.17 

.100 

6.096-2 

ELEV 

4.106-5 

4.116-2 

5.716-2 

DELEV 

6.17 

3.816-2 

5.496-2 

UM 

-2.576-3 

-5.386-2 

—6.746—2 

RHL4 

.131 

6.386-2 

5.976-2 

D6VL 

-8.236-3 

-3.876-2 

-5.396-2 
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Table  7  Regression  Results  for  Case  #27  (Box  30,  Fall,  21^  9-Hour 
Forecast)  Showing  Top  15  Predictors,  Regression  Weights 
(Raw  and  Standardized),  and  Partial  Correlations 


VARIABLE  NAME 

B 

BSTD 

PARTIAL 

CORRELATION 

ACCL9 

.229 

.176 

.133 

RHM 

.179 

9.54e-2 

7.27e-2 

NEPH9 

.150 

.137 

9.39e-2 

NEPH21 

6.50e-2 

6.51  e-2 

6.70e-2 

ZENITH 

-.201 

-6.46e-2 

— 8.68e— 2 

UM 

-4.04e-3 

-8.71  e-2 

-.107 

ACCM9 

.141 

.108 

7.61  e-2 

NEPH18 

5.91  e-2 

5.97e-2 

5.53e-2 

D9VTL 

1.63e3 

4.59e-2 

5.82e-2 

DELEV 

6.30 

3.89e-2 

5.39e-2 

ACCH9 

8.73e-2 

6.67e-2 

6.00e-2 

ELEV 

5.01  e-5 

3.91  e-2 

5.24e-2 

RH90 

5.62e-2 

5.88e-2 

5.99e-2 

09RHM 

.145 

5.10e-2 

5.83e-2 

NEPH12 

6.81  e-2 

6.44e-2 

4.65e-2 
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Table  8  Regression  Results  for  Case  #42  (Box  44,  Winter,  OZ,  3-Hour 
Forecast)  Showing  Top  15  Predictors,  Regression  Weights 
(Raw  and  Standardized),  and  Partial  Correlations 


VARIABLE  NAME 

B 

BSTD 

PARTIAL 

CORRELATION 

ACCL3 

.288 

.238 

.167 

NEPH3 

.341 

.336 

.330 

ACCH3 

.101 

8.46e-2 

6.67e-2 

RHM2 

6.47e-2 

3.81  e-2 

4.24e-2 

ACCM3 

.111 

9.26e-2 

5.65e-2 

TM 

2.90e-2 

4.99e-2 

7.22e-2 

VTL 

6.54e+2 

3.40e-2 

5.04e-2 

NEPH18 

2.87e-2 

3.03e-2 

3.91  e-2 

TAH 

1.37e+2 

3.97e-2 

5.72e-2 

RH70 

3.38e-2 

3.98e-2 

4.93e-2 

D3P 

-1.00e-2 

-3.34fr-2 

-4.80e-2 

NEPH24 

3.27e-2 

3.31  e-2 

4.37e-2 

D3RHL 

.187 

2.90e-2 

4.19e-2 

TYPE 

1.43e-2 

1.92e-2 

3.00e-2 

LON 

-8.00e-2 

-1.84e-2 

-2.80e-2 
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'Dable  9  Regression  Results  for  Case  #52  (Box  61,  Summer,  15Z, 
6-Hour  Forecast)  Showing  Top  IS  Predictors,  Regression  Weights 
(Raw  and  Standardized),  and  Partial  Correlations 


VARIABLE  NAME 

B 

BSTD 

PARTIAL 

CORRELATION 

ACCH6 

.149 

.133 

9.25e-2 

LAT 

-2.51  e-2 

-.233 

-.158 

NEPH24 

9.08e-2 

8.97e-2 

9.19e-2 

ACCM6 

.121 

.107 

7.22e-2 

WM 

9.75e3 

.159 

.118 

VTL 

2.27e3 

6.66e-2 

6.57e-2 

LNRHM 

4.61  e-2 

4.19e-2 

3.64e-2 

D6WH 

-9.95e3 

-6.44e-2 

-6.50e-2 

D6UL 

1.17e-2 

5.06e-2 

5.05e-2 

D6P 

-1.49e-4 

-3.85e-2 

-4.03e-2 

RHAB 

9.84e-2 

3.70e-2 

3.44e-2 

NEPH6 

3.80e-2 

3.81  e-2 

3.22e-2 

DELEV 

1.75 

2.86e-2 

3.17e-2 

VM 

6.45e-3 

3.20e-2 

3.21  e-2 

UL 

3.30e-3 

3.70e-2 

2.61  e-2 
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Partial  Correlations 


One  way  to  analyze  the  results  of  regression  modeling  is  to  look  at  partial  correla¬ 
tions.  Partial  correlations  are  related  to  the  contribution  of  each  individual  independent 
variable  to  the  prediction  of  the  dependent  variable  (cloud  cover)  after  predictability 
derived  from  previously  selected  independent  variables  has  been  accounted  for.  The 
squared  partial  correlation  is  a  measure  of  that  part  of  the  residual  variance  accounted 
for  by  each  variable. 

The  partial  correlations  for  the  five  sample  cases  are  hsted  along  with  the  regres¬ 
sion  weights  in  Tables  5  through  9.  Note  that  in  the  tables,  variables  are  ordered  according 
to  their  order  of  selection  from  the  original  set  of  approximately  80  predictors.  \\^thin  the 
fifteen  predictors  in  each  table,  selection  order  would  have  been  from  the  largest  to  small¬ 
est  (magnitude)  particil  correlation. 

Reduction  in  Variance 

Another  way  to  study  the  results  of  regression  modeling  is  to  look  at  the  reduction 
in  variance  due  to  the  individual  predictors  and  due  to  the  complete  regression  model.  We 
can  evaluate  the  MOS  approach  to  cloud  amount  prediction  by  comparing  the  reduction 
in  variance,  R^,  from  the  MOS  equations  to  that  obtained  using  only  the  persisted  cloud 
cover  amount.  That  comparison  (organized  into  three  tables  by  box  number),  for  the  71 
cases  we  modeled,  is  shown  in  Tables  10a  through  10c. 
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Table  10a  Comparison  of  the  Reduction  in  Variance  (R^)  Due  to 
MOS  and  Persistence  for  Box  30  Cases — Values  Marked  With  *‘L’’ 
and  Correspond  to  Land  and  Water  Cases,  Respectively 


VALID 

TIME 

(Z) 

FORECAST 

LENGTH 

(Hours) 

SPRING 

SUMMER 

FALL 

WINTER 

MOS 

Persist 

MOS 

Persist 

MOS 

Persist 

MOS 

Persist 

3 

3 

.5370 

.3569 

3 

6 

.4909 

.2761 

3 

9 

.4487 

.2012 

9 

3 

.5831 

.4638 

.6222 

.4946 

.6430 

.5336 

.5672 

.4381 

9 

6 

.4868 

.2982 

.5434 

.3494 

.5538 

.3649 

.4831 

.3052 

9 

6 

.5005 

.2982 

.5614 

.3494 

.5563 

.3649 

.4780 

.3052 

9 

9 

.4468 

.1761 

.5221 

.2489 

.5043 

.2534 

.4199 

.2086 

21 

3 

.5230 

.4048 

.5510 

.4572 

.6096 

.5010 

.5385 

.4295 

21 

6 

.4132 

.2293 

.4448 

.2776 

.5283 

.3484 

.4300 

.2651 

21 

6 

.4013 

.2293 

.4297 

.2776 

.5194 

.3484 

.4212 

.2651 

21 

9 

.3661 

.1627 

.4002 

.2232 

.4921 

.2938 

.4063 

.2092 

21 

3 

.5378'- 

.4315*- 

21 

6 

.42861- 

.2587'- 

21 

9 

.3916'- 

.1972'- 

21 

3 

.5311'"' 

.3831'''' 

21 

6 

.3829'"' 

.1519'« 

21 

9 

.3808'''' 

.1347'''' 

Table  10b  Comparison  of  the  Reduction  in  Variance  (R^)  Due  to 
MOS  and  Persistence  for  Box  44  Cases 


VALID  FORECAST 
TIME  LENGTH 

(Z)  (Hours) 


3 


6 


9 


3 


6 


9 


SPRING 


MOS  Persist 


SUMMER 


MOS  Persist 


FALL 


MOS  Persist 


WINTER 


.5062 

.3519 

.4181 

.1910 

.3727 

.1109 

.5378 

.3513 

.4250 

.1760 

.3938 

.1288 

MOS 

Persist 

.6124 

.4784 

.5340 

.2300 

.4636 

.1890 

.5258 

.3785 

.4448 

.2192 

.4247 

.1954 

Table  10c  Comparison  of  the  Reduction  in  Variance  (R^)  Due  to 
MOS  and  Persistence  for  Box  61  Cases  —  Values  Marked  With  **1." 
and  Correspond  to  Land  and  Water  Cases,  Respectively 


VALID 

TIME 

(Z) 

FORECAST 

LENGTH 

(Hours) 

SPRING 

FALL 

WINTER 

MOS 

Persist 

MOS 

Persist 

MOS 

Persist 

MOS 

Persist 

3 

3 

.4168 

.2549 

.5594 

.4186 

3 

6 

.4378 

.0688 

.4055 

.1761 

3 

9 

.1677 

.0518 

.3843 

.1426 

3 

3 

.5324*- 

.4018'- 

3 

6 

.3525^ 

3 

9 

.3377L 

.1628'- 

3 

3 

.5693'^ 

3 

6 

.5030'^ 

.1472'^ 

3 

9 

.5263^ 

.1050'*' 

15 

3 

.3387 

.2106 

.3979 

.2684 

15 

6 

.1940 

.0721 

.2930 

.1693 

15 

9 

.1830 

.0625 

.2884 

.1619 

From  the  entire  71  cases  we  were  able  to  draw  a  few  conclusions. 


•  The  reduction  in  variance  due  to  MOS  or  persistence  decreases  with  increasing 
forecast  length.  That  is,  for  a  3-hour  forecast  is  greater  than  for  a  6-hour 
forecast  which  is  in  turn  greater  than  R^  for  a  9-hour  forecast  for  a  given  box, 
season,  and  time. 

•  The  percentage  increase  in  the  reduction  of  variance  using  the  MOS  approach 
over  that  using  persistence  varied  substantially  from  one  region  and  season  to 
the  next.  The  greatest  increases  were  found  in  box  61  summer  3Z  and  winter  3Z 
(most  notable  was  the  box  61,  winter,  3Z,  9-hour  forecast  for  water  grid  points). 

•  R2  due  to  persistence  was  higher  in  the  winter  than  the  summer  in  boxes  44  and 
61.  Its  highest  values  were  in  box  30,  fall.  Lowest  values  in  box  61,  summer. 

•  The  most  variability  in  R^  values  was  found  in  the  box  61  cases. 

In  all  of  these  cases,  the  MOS  approach  improved  upon  persistence  because  it  used 
more  information.  Each  additional  variable  that  is  used  in  the  regression  equation  adds 
to  the  overall  reduction  in  variance,  though  each  successive  predictor  adds  less  and  less 
to  the  overall  reduction,  in  general. 

In  addition  to  stud3dng  the  overall  reduction  in  variance  due  to  all  fifteen  predictors 
in  each  case  we  also  studied  the  incremental  changes  in  R^  due  to  individual  predictors. 
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Table  11^  shows  the  cumulative  reduction  in  variance  for  two  sample  cases  (numbers  42 
and  52)  as  successive  variables  were  added  to  the  regression  equations.  From  the  results 
presented  in  Table  11  it  is  clear  that  the  bulk  of  the  variance  reduction  comes  from  the  first 
few  predictors.  In  terms  of  variance  reduction,  our  decision  to  retain  fifteen  predictors  in 
the  models  was  very  conservative. 


Table  1 1  Incremental  Changes  to  for  Cases  #42  (Box  44,  Winter, 

OZ,  3-Hour  Forecast)  and  #52  (Box  61,  Summer,  15Z,  6-Hour  Forecast) 


Case#42  Box  44 

Case#52  Box  61 

Variable 

Cumulative 

Variable 

Cumulative 

ACCL3 

.5187788 

ACCH6 

.1166509 

NEPH3 

.5781473 

LAT 

.1387272 

ACCH3 

.5920656 

NEPH24 

.1539613 

RHM2 

.5962372 

ACCM6 

.1637530 

TM 

.5979326 

WM 

.1706838 

TAH 

.5996365 

VTL 

.1780713 

NEPH24 

.6011909 

LNRHM 

.1818255 

RH70 

.6025454 

D6WH 

.1851151 

ACCM3 

.6037927 

D6UL 

.1853130 

D3P 

.6048484 

D6P 

.1890684 

VTL 

.6057851 

RHAB 

.1901213 

D3RHL 

.6064843 

NEPH6 

.1910302 

NEPH18 

.6071277 

DELEV 

.1918519 

TYPE 

.6073927 

VM 

.1926227 

LON 

.6076998 

UL 

.1931708 

Residual  Analysis 

A  final  step  in  our  analysis  was  to  examine  the  residuals.  Linear  regression  as¬ 
sumes  linear  relationships  between  each  of  the  independent  variables  and  the  predictand. 


^  The  values  included  in  Tables  10a  through  10c  correspond  to  the  reduction  in  variance 
due  to  the  strongest  fifteen  predictors  of  the  pool  of  approximately  80  predictors  used 
during  the  forward  regression  step.  The  values  reported  in  Table  11  were  determined 
during  the  regression  modeling  step  when  only  the  strongest  fifteen  predictors  were  used. 
Slight  discrepancies  between  the  numbers  are  due  to  the  influence  of  the  remaining 
predictors  not  used  in  the  final  regression  modeling  step. 


73 


Linear  regression  also  assumes  that  the  residual  population  density  is  unimodal  and  sym* 
metric  around  zero.  By  analyzing  a  “normal  plot  of  the  residuals”  (one  of  the  options  avail¬ 
able  through  STATISTIC  A™ )  we  were  able  to  verify  whether  our  model  residuals  were 
normal.  Testing  for  normality  is  a  stronger  test  than  required  by  the  assumptions  of  linear 
regression,  but  it  can  point  to  potential  problems  in  the  data. 

A  normal  probability  plot  of  the  residuals  is  generated  by  plotting  the  actual  residu¬ 
als  on  the  X  axis  and  the  z-values  (standardized  assuming  a  normal  distribution)  of  the 
rank-ordered  residuals  on  the  y  axis.  If  the  residuals  follow  a  normal  distribution,  the  data 
lie  along  a  straight  line  (the  “normal”  line).  Gross  deviations  hum  normahty  show  up  as 
deviations  from  the  normal  line.  In  certain  cases  some  non-linearities  can  be  removed  by 
transforming  the  data.  Normal  plots  of  the  residuals  for  the  five  sample  cases  mentioned 
above  are  included  in  Figures  36  through  40. 

In  four  of  the  cases  that  we  looked  at,  the  data  did  not  show  any  large  deviations 
from  the  ideal  normal  line.  Significant  deviations  were  found  in  box  61,  shown  in  Fig¬ 
ure  40.  The  most  negative  residual  points  that  fall  well  below  the  normal  line  indicate  a 
heavy-tailed  distribution.  The  flattening  of  the  curve  around  residual  values  between  -0.5 
and  0  indicates  a  tendency  toward  a  multi-modal  density  function.  A  histogram  of  the  re¬ 
siduals  for  the  box  61  case  (see  Figure  41)  clearly  supports  both  of  those  indications.  In 
contrast.  Figure  42  shows  a  histogram  of  the  model  residuals  from  the  case  presented  in 
Figure  38  that  fell  nicely  along  the  normal  line.  In  this  case,  the  residual  distribution 
clearly  satisfies  the  zero  mean  and  unimodal  requirements  assumed  in  linear  regression. 

3.4.2  Predicting  Probabilities  of  Cloud  Cover  Categories 

An  introduction  to  the  process  of  predicting  probability  of  cloud  cover  categories 
was  given  in  Section  3. 1.  We  used  the  actual  error  distributions  resulting  from  the  regres¬ 
sion  modeUng  to  determine  those  probabilities.  First  we  divided  the  predicted  cloud 
amounts  into  11  bins:  0-4%,  5-14%,  . . .  ,  85-94%,  95-100%.  Then  we  generated  histo¬ 
grams  of  the  errors  (binned  in  0.1  intervals)  by  forecast  cloud  category.  In  general,  these 
residuals  have  values  between  -1  and  +1.  For  any  given  cloud  amount  estimate,  y,  howev¬ 
er,  the  residual  values  can  range  from  a  low  of  -  y  to  a  high  of  1-y.  Due  to  the  bounded 
nature  of  cloud  cover  observations,  the  error  histograms  frequently  followed  U-,  J-,  or  L- 
shaped  distributions.  Figure  43  shows  the  residual  distributions  for  each  of  the  forecast 
categories  for  case  #27  (box  30/Fall/21Z/9-hour  forecast).  A  histogram  of  the  entire  residu¬ 
al  data  set  for  this  case  was  shown  in  Figure  42. 
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As  shown  in  Section  3.1,  the  probability  that  the  observed  cloud  cover  will  be  less 
than  or  equal  to  a  threshold  value,  c,  is  given  by 


PriY  ^  c  I  Yl  *  PrIE  s  c-Y  I  Y| 


(3-17) 


where  Y  is  the  observed  cloud  cover,  Y  is  the  estimated  cloud  cover  using  the  MOS  equa¬ 
tions,  and  E  is  the  error.  For  example,  to  find  the  probability  that  the  observed  cloud  cover 
will  be  ^  50%  given  that  the  predicted  cloud  amount  is  80%,  we  determined  the  ctimula- 
tive  density  (i.e.,  probability)  below  -0.3  (i.e.,  0.5-0.8)  in  the  error  histogram  for  cloud  cov¬ 
er  forecasts  in  the  range  75-84%.  The  probability  of  a  single  cloud  cover  category  (i.e.,  that 
the  observed  cloud  cover  is  between  two  thresholds)  is  simply  the  difference  between  the 
probabilities  computed  for  each  of  the  thresholds  individually. 

These  calculations  do  not  make  any  assumptions  about  the  behavior  of  the  errors. 
There  is  some  approximation  involved,  however,  due  to  the  finite  widths  of  the  bins  used 
to  categorize  cloud  cover  forecasts  and  the  errors  themselves. 

We  generated  MOS-based  probabilities  of  cloud  cover  categories  at  10%  intervals 
for  all  7 1  cases  listed  in  Table  3.  Figures  44  through  48  show  the  results  for  the  five  sample 
cases  discussed  previously. 


3.4.3  Modeling  with  the  KEEP  Methodology 

In  Section  3. 1  we  introduced  the  regression  estimation  of  event  probabilities  meth¬ 
odology  to  estimate  cloud  cover  categories.  In  this  section  we  discuss  how  we  implemented 
the  KEEP  approach,  show  results  of  a  sample  case  (box  44/winter/3Z/3-hour  forecast),  and 
discuss  why  we  decided  to  implement  an  alternate  methodology  to  predict  both  total  cloud 
cover  amount  and  probability  of  cloud  cover  categories. 

First,  we  transformed  the  response  variable,  total  cloud  cover  (cc),  to  6  categorical 
response  variables  defined  in  consultation  with  Phillips  Laboratoiy  personnel  as 


Category 

1: 

cc  =  0% 

Category 

2; 

0  <  cc  ^  25% 

Category 

3: 

25  <  cc  ^  50% 

Category 

4; 

50  <  cc  ^  75% 

Category 

5: 

75  <  cc  <  100% 

Category 

6: 

cc  =  100%. 

The  result  was 

a  i 

binary  response  \ 

1  one. 
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The  KEEP  methodology  uses  a  common  predictor  set  to  build  a  regression  equation 
for  each  category;  this  ensures  that  the  probabilities  sum  to  one.  Therefore,  the  next  step 
was  to  determine  a  set  of  predictors  common  to  all  of  our  six  categories.  We  found  that  this 
step  was  somewhat  subjective.  We  used  the  forward  stepwise  regression  module  to  select 
the  strongest  predictors  for  each  category  individually.  Then  we  made  a  histogram  of  the 
strongest  predictors  across  all  six  categories  and  selected  those  variables  that  occurred 
most  frequently.  When  two  or  more  variables  occurred  with  the  same  frequency,  we  se¬ 
lected  that  one  which  was  the  “stronger”  predictor  (i.e.,  occurred  closer  to  the  top  of  the  list 
of  predictors).  Table  12  contains  a  list  of  the  fifteen  strongest  predictors  for  each  of  the  six 
categories  as  selected  using  forward  regression.  Notice  that  in  one  of  the  categories  only 
8  predictors  were  selected  (i.e.,  only  those  eight  were  found  to  be  significant  using  the  F 
test  provided  in  the  STATISTICA"*  model  software).  The  last  column  contains  those  vari¬ 
ables  we  chose  for  the  common  predictor  set  in  the  sample  case. 

Next,  we  used  KEEP  to  develop  six  models,  one  for  each  of  the  six  categories,  using  the 
fifteen  common  predictors.  The  regression  coefficients  for  each  category  in  the  sample  case 
are  listed  in  Table  13  (we  list  both  raw  and  standardized  weights  as  in  the  previous  section). 


Table  12  List  of  Strongest  Predictors  Selected  for  Each  of  the  Six 
Cloud  Cover  Categories  and  the  Conunon  Predictor  Set  Used  in  the 
KEEP  Case  (Box  44,  Winter,  3Z,  3-Hour  Forecast) 


CATEGORY 

1 

CATEGORY 

2 

CATEGORY 

3 

CATEGORY 

4 

CATEGORY 

5 

CATEGORY 

6 

COMMON 
SET  OF 
PREDICTORS 

1 

ACCL3 

TDIF1 

LAT 

NEPH3 

ACCL3 

NEPH3 

NEPH3 

2 

NEPH3 

VTL 

DAY 

ELEV 

ELEV 

ACCL3 

ACCL3 

3 

RHM 

RH50 

RHL2 

MCC3 

LAT 

RHM2 

LON 

4 

NEPH21 

RHL2 

NEPH3 

UL 

LON 

ZENITH 

RHL2 

5 

TL 

ZENITH 

TYPE 

LNRHM 

DAY 

NEPH21 

ELEV 

6 

ACCH3 

WH 

UL 

TDIF1 

SFCP 

TN 

TYPE 

7 

TYPE 

LNRHM 

D3VL 

TYPE 

DELEV 

WH 

SFCP 

8 

ELEV 

RHM2 

SFCP 

DELEV 

TAH 

TAH 

RH50 

9 

RH50 

LON 

SFCP 

D3WL 

RHM 

TDIF1 

10 

D3VTH 

WRHM 

LON 

NEPH3 

ACCM3 

ZENITH 

11 

TAM 

RH90 

RH50 

NEPH6 

NEPH9 

NEPH21 

12 

LON 

DVM 

RHL2 

D3RHL 

RHL2 

UL 

13 

TDIF1 

TM 

ACCL3 

D3VTL 

LAT 

RHM2 

14 

D3WL 

NEPH18 

RHL 

TAL 

NEPH15 

RHM 

15 

D3RHM 

VTL 

RHAB 

TM 

TAM 

LNRHM 

Table  13  Results  for  REEP  Case  (Box  44,  Winter,  3Z,  3-Hour  Forecast)  Showing  the  Top  15 
Predictors  Common  to  All  Categories  and  Corresponding  Egression  Weights  (Raw  and  Standardized) 
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The  equations  predict  probabilities  of  cloud  cover  categories,  but  unlike  true  proba¬ 
bilities,  these  individual  values  are  not  required  to  be  between  zero  and  one.  For  example, 
the  KEEP  model  produces  the  following  probability  estimates  at  a  single  grid  point  on  the 
first  day  of  the  sample  case  (Dec.  1,  box  44): 

Category  1:  -0.120664 

Category  2:  0. 1000 10 

Category  3:  0.013345 

Category  4:  0.221180 

Category  5:  0.054954 

Category  6:  0.731176. 

The  probabilities  sum  to  one  as  required,  but  are  not  bounded  by  0  and  1  as  a  result 
of  using  a  non-bounded  linear  regression  model.  The  predicted  total  cloud  amount  can  be 
determined  from  the  probabiUties  using  a  prescribed  selection  criterion  such  as  the  most 
probable  category  or  the  median  category.  In  the  example  above,  either  of  these  criteria 
would  result  in  a  cloud  cover  estimate  of  category  6;  100%. 

Besides  the  two  disadvantages  described  in  Section  3.1, 

1.  high  computational  cost  to  determine  a  common  predictor  set  for  all  six  categorical 
response  variables 

2.  probabilities  outside  the  range  [0,1], 

our  trial  test  showed  that  the  data  would  not  support  predicting  middle  (partly-cloudy) 
categories.  This  is  evident  when  looking  at  the  reduction  in  variance,  R^,  for  the  six  cloud 
cover  categories  shown  in  Table  14.  The  reduction  in  variance  using  the  KEEP  approach 
is  relatively  high  for  the  two  outer  (clear  and  overcast)  categories  and  very  low  for  inter¬ 
mediate  categories.  There  are  very  few  data  in  these  categories  and  thus  little  to  support 
model  development. 

Table  14  Results  of  KEEP  Case  Showing  Reduction  in  Variance  (R^) 
for  the  Six  Categorical  Cloud  Cover  Regression  Models.  Note  the 
Low  Values  for  the  Four  Middle  Categories 
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3.5  ESTIMATING  CURRENT  AND  FUTURE  MODEL  PERFORMANCE 

How  does  the  MOS  approach  compare  to  current  technology?  How  will  the  equa¬ 
tions  developed  under  this  effort  perform  in  the  future  with  independent  weather  data? 
This  section  seeks  to  answer  both  of  those  questions. 

First,  we  compare  the  MOS  approach  to  persistence.  This  is  not  an  arbitrary  refer¬ 
ence  point.  Persistence  is,  in  general,  a  good  predictor  of  cloud  cover  and  is  relied  on  heavi¬ 
ly  in  certain  parts  of  the  world  as  a  forecasting  tool.  For  example,  the  TRONEW  and 
HRCP  models  both  use  a  type  of  persistence  to  model  cloud  cover.  Second,  we  estimate 
future  performance  using  the  jackknife  technique  which  subsamples  the  model  develop¬ 
ment  data  set,  develops  independent  models  for  each  subset,  then  assesses  performance 
variability  when  each  model  is  applied  to  data  not  used  in  its  development.  Finally,  we 
validate  the  cloud  cover  category  forecasts  by  developing  a  model  based  on  a  subset  of  the 
development  data  set  and  comparing  those  model  results  to  the  results  computed  by  ap¬ 
plying  the  same  model  to  the  remainder  of  the  data. 

3.5.1  Comparison  of  MOS  Approach  to  Persistence 

To  answer  the  question  “How  does  the  MOS  approach  compare  to  persistence?”  one 
can  use  any  one  of  a  number  of  statistical  performance  statistics.  These  include:  the  F  sta¬ 
tistic  for  significance  of  the  regression,  the  residual  mean  square,  the  coefficient  of  deter¬ 
mination,  the  total  squared  error,  etc.  To  compare  the  two  forecast  approaches  then,  one 
compares  performance  statistics.  We  chose  to  use  performance  statistics  more  familiar  to 
the  weather/forecasting  community  than  those  listed  above  to  facilitate  comparison  with 
other  forecast  models. 

We  selected,  in  conjunction  with  Phillips  Laboratory  personnel,  three  performance 
statistics:  the  Brier  score,  the  20/20  score,  aiid  sharpness.  Sharpness  is  not  strictly  a  per¬ 
formance  statistic  (i.e.,  it  is  not  a  measure  of  the  bivariate  forecast-observation  distribu¬ 
tion),  rather  it  is  a  measure  of  a  single  cloud  cover  distribution  (e.g.,  of  forecast  data), 
which  can  be  compared  to  other  distributions  (e.g.,  of  observed  data).  That  is,  sharpness 
is  a  measure  of  the  forecast  and  observed  marginal  distributions. 

We  begin  our  discussion  with  a  description  of  the  21x21  bin  contingency  tables 
used  to  contain  the  forecast  and  observed  data  from  which  we  compute  the  performance 
measures. 
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The  Contingency  Table 

A  contingency  table  is  simply  a  two-dimensional  histogram  with  observed  cloud 
cover  along  one  axis  and  forecast  cloud  cover  (generated  by  MOS  or  persistence)  along  the 
other  axis.  We  used  a  standard  binning  convention  that  assigns  cloud  cover  (cc)  to  21  in¬ 
crements  as  follows: 

bin  1:  cc  <  2.5% 

bin  2:  2.5%  ^  cc  <  7.5% 

bin  3:  7.5%  <  cc  <  12.5% 


bin  19:  87.5%  s  cc  <  92.5% 

bin  20:  92.5%  ^  cc  <  97.5% 

bin  21:  cc  >  97.5%. 

The  resulting  contingency  table  has  21x21  elements  (or  cells).  A  schematic  of  a 
21x21  contingency  table  is  provided  in  Figure  49.  For  each  of  the  71  cases  we  analyzed, 
we  built  contingency  tables  containing  one  season  of  data  (approximately  90  days)  for  225 
eighth-mesh  grid  points  ( 15  x  15  points),  a  maximum  of  approximately  20,250  counts  (box 
61  cases  had  fewer  points  because  a  portion  of  the  box  contains  “off  world”  points,  see  Fig¬ 
ure  3).  The  total  number  of  counts  in  the  table  is  denoted  by  N.  In  most  cases,  N  was  much 
less  than  the  maximum  after  accounting  for  missing  and  “old”  cloud  cover  data  points.  The 
horizontal  (i)  axis  of  the  table  in  Figure  49  corresponds  to  the  observed  cloud  cover.  The 
vertical  (j)  axis  of  the  table  corresponds  to  the  forecast  cloud  cover.  For  any  (i,  j)  cell,  the 
number  of  counts  in  the  cell  is  denoted  by  ny. 

The  table  is  filled  one  count  at  a  time  by  determining  the  cell  that  corresponds  to 
the  forecast  value  along  the  row  axis  and  observed  cloud  cover  along  the  column  axis.  If 
the  forecast  is  perfect  at  every  point  over  the  entire  season,  all  of  the  entries  in  the  contin¬ 
gency  table  will  lie  along  the  diagonal  where  forecast  cloud  cover  is  equal  to  observed 
cloud  cover. 

The  Brier  score,  20/20  score,  and  sharpness  are  measures  of  how  the  counts  are  dis¬ 
tributed  throughout  the  contingency  table.  They  are  described  below. 
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The  Brier  Score 


The  Brier  score  (Ref.  3),  measures  the  mean  squared  difference  between  the  fore¬ 
cast  (F)  and  observed  (O)  cloud  cover  amounts  (each  ranges  from  0.0  to  1.0)  from  the 
21x21  contingency  table.  The  Brier  score,  Ssnen  is  computed  by  summing  over  all  data 
points  in  a  season  as  follows: 

21  21 

Sfirier  =  H  2  X  -  Ol'"' 

j=l  i=l 


This  score  is  particularly  sensitive  to  counts  in  the  extreme  off-diagonal  cells  due  to  the 
quadratic  term.  The  Brier  score  ranges  from  a  perfect  (minimum)  value  of  zero,  corre¬ 
sponding  to  all  data  points  lying  along  the  diagonal,  to  a  maximum  of  one,  corresponding 
to  all  data  points  lying  in  the  extreme  bins. 

We  computed  Brier  scores  for  the  MOS-based  forecast  and  the  persistence  forecast 
for  each  of  the  model  cases.  Since  the  Brier  score  measures  mean  squared  error  and  linear 
least  squares  regression  minimizes  squared  error,  we  expected  that  the  MOS-based  fore¬ 
casts  would  perform  well  by  this  metric.  Indeed,  MOS  outperformed  persistence  in  all  of  the 
71  cases.  The  scores  are  listed  in  Table  15.  Figures  50  through  52  show  bar  graphs  of  the 
Brier  scores  for  each  of  the  RTNEPH  boxes  that  we  studied:  30,  44,  and  61,  respectively. 

The  cyclic  behavior  evident  in  Figures  50  through  52  is  caused  by  ordering  the  mod¬ 
el  cases  by  forecast  length  (refer  to  Table  3).  The  cases  are  listed  in  order  of  3-,  6-  (some 
cases  have  two  6-hour  forecasts),  and  9-hour  forecasts.  By  looking  at  those  cases  in  Box 
30  for  which  we  generated  two  6-hour  forecasts  (one  each  for  the  OZ  and  12Z  initialized 
GSM  data),  we  see  that  there  is  no  evidence  to  suggest  that  one  initialization  time  is  in¬ 
herently  better  than  another.  This  implies  that  the  MOS  equations  compensated  for  any 
model  spin-up  problems  that  may  have  existed  in  the  GSM. 

If  we  explicitly  average  all  the  forecasts  of  a  given  length  from  one  box  and  compare 
the  MOS-based  and  persistence-based  average  scores  for  differing  forecast  lengths  (see 
Figure  53),  it  is  clear  that  Brier  scores  increase  with  increasing  forecast  length.  That  is, 
Sfirier  for  a  9-hour  forecast  is  greater  than  Sfiner  for  a  6-hour  forecast  which  is  greater  than 
SBrier  for  a  3-hour  forecast  (for  both  persistence  and  MOS).  This  implies  the  obvious,  that 
shorter  forecasts  are  more  skillful  (i.e.,  have  less  error)  than  longer  forecasts. 
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Table  15  Brier  Scores,  20/20  Scores,  and  Sharpness  for  the  7 1  Model 

Cases  Introduced  in  Table  3 


c 

ASE 

N. 

AME 

30/SP/3/3/12 


2  30/SP/3/6/12 


3  30/SP/3/9/12 


30/SP/9/3/0 


30/SP/9/6/0 


30/SP/9/6/12 


30/SP/9/9/12 


30/SP/21/3/12 


9  30/SP/21/6/12 


10  30/SP/21/«5/0 


30/SP/2 1/9/0 


30/SU/9/3/0 


3(VSU/9/6/0 


30/SU/9/6/12 


30/SU/9/9/12 


30/SU/21/3/12 


30/SU/21/6/12 


3C/SU/2 1/6/0 


30/SU/2 1/9/0 


30/FA/9/3/0 


30/FA/9/6/0 


22  30/FA«/6/12 


23  30/FA/9/9/12 


30/FA/21/3/12 


30/FA/21/6/12 


26  30/FA/21/6/0 


30/FA/21/9/0 


30/WI/9/3/0 


30/WI/9/6/0 


30/W1/9/6/12 


30/WI/9/9/12 


30/WI/21/3/12 


30/WI/21/6/12 


30/WI/21/6/0 


30/WI/21/9/0 


44/SU/0/3/12 


44/SU/0/6/12 


44/SU/0/9/12 


44/SU/1 5/3/0 


44/SU/1 5/6/0 


0.0942 


0.1017 


0.1092 


0.0679 


0.0833 


0.0839 


0.0932 


0.0909 


0.1112 


0.1143 


0.1216 


0.0589 


0.0714 


0.0712 


0.0777 


0.0723 


.0915 


0.0939 


0.1018 


0.0620 


0.0777 


0.0785 


0.0884 


0.0745 


0.0918 


0.0937 


0.0999 


0.0741 


0.0900 


0.0929 


0.104^ 


0.0921 


0.1138 


0.1157 


0.1216 


0.0631 


0.0754 


0.0827 


0.0840 


0.1055 


PERSISTED 


0.1684 


0.1865 


0.2011 


0.1115 


0.1566 


0.1566 


0.2200 


0.1311 


0.1835 


0.1838 


0.2090 


0.1001 


0.1484 


0.1464 


0.1801 


0.1069 


0.1588 


0.1583 


0.1808 


0.0980 


0.1401 


0.1409 


0.1901 


0.1075 


0.1438 


0.1439 


0.1591 


0.1163 


0.1568 


0.1565 


0.2038 


0.1341 


0.1758 


0.1765 


0.1949 


0.1021 


0.1505 


0.2083 


0.1500 


0.2065 


20/20 

MOS 

20/20 

PERSiSTEO 

SHARPNESS 

MOS 

SHARPNESS 

PERSISTEO 

0.5651 

0.6238 

0.4877 

0.8007 

0.5365 

0.5911 

0.4632 

0.7925 

04821 

0.5346 

0.4084 

0.7111 

0.6767 

0.6644 

0.5499 

0.7977 

0.5956 

0.5852 

0.4426 

0.8222 

0.5846 

0.4531 

0.8201 

0.5285 

0.5016 

0.3653 

0.8035 

0.5867 

0.6201 

0.4988 

0.7030 

0.4677 

0.5149 

0.3627 

0.6330 

0.4635 

0.5129 

0.3602 

0.6311 

0.4286 

0.4828 

0.3221 

0.6301 

0.7129 

0.6733 

0.5794 

0.8062 

0.6354 

0.5859 

0.4825 

0.8246 

0.6418 

0.5887 

0.4925 

0.8251 

0.5962 

0.5374 

0.4280 

07992 

0.6549 

0.6420 

0.5637 

0.6908 

0.5550 

0.5263 

0.4519 

0.6071 

0.5431 

0.5293 

0.4424 

0.6086 

0.5098 

0.4932 

0.4061 

0.5953 

0.7131 

0.7098 

0.6100 

0.8104 

0.6299 

0.6373 

0.5305 

0.8480 

0.6302 

0.6335 

0.5308 

0.8474 

0.5742 

0.4556 

0.8421 

0.6739 

0.6999 

0.7864 

0.5881 

0.6037 

0.5144 

0.7051 

0.5797 

0.6061 

0.5087 

0.7082 

0.5530 

0.5804 

0.4875 

0.6984 

0.6698 

0.7055 

0.5972 

0.8327 

0.5826 

0.6429 

0.4935 

0.8592 

0.5813 

0.6457 

0.5044 

0.8632 

0.5(»7 

0.5920 

0.4080 

0.8562 

0.6045 

0.6815 

0.5504 

0.8152 

0.4970 

0.5946 

0.4373 

0.7560 

0.4932 

0.5956 

0.4309 

0.7571 

0.4590 

0.5682 

0.3927 

0.7555 

0.6438 

0.5749 

0.3588 

0.3616 

0.5596 

0.4525 

0.2790 

0.4401 

0.5193 

0.3847 

0.2367 

0.6783 

0.6143 

0.5976 

0.5315 

0.7466 

0.4964 

0.4750 

0.3971 

0.6588 

SHARPNESS 


0.8222 


0.8251 


0.8294 


0.7076 


0.7037 


0.7031 


0.7026 


0.7896 


0.7929 


0.7956 


0.7991 


0.7028 


0.6913 


0.6927 


06855 


0.7835 


0.7798 


0.7809 


0.7792 


0.7560 


0.7583 


0.7569 


0.7586 


0.8200 


0.8260 


0.8269 


0.8327 


0.7996 


0.7982 


0.8011 


0.8016 


0.8384 


0.8450 


0.8451 


0.8433 


0.3917 


0.3990 


0.4077 


0.7230 


0.7231 
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Table  15  Brier  Scores,  20/20  Scores,  and  Sharpness  for  the  7 1  Model 
Cases  Introduced  in  Ihble  3  (Continued) 


F 

CASE 

NAME 

N 

BRIER 

MOS 

BRIER 

PERSISTED 

20/20 

MOS 

20/20 

PERSISTED 

SHARPNESS 

PERSISTED 

SHARPNESS 

OBSERVED 

IB 

44/SU/15«/0 

12609 

0.1118 

0.2325 

0.4644 

0.4601 

0.3620 

0.7066 

0.7274 

44/WI/(V3/12 

12450 

0.0675 

0.1044 

0.6746 

0.6225 

0.5369 

0.5815 

0.6223 

la 

44/WI/0/6/12 

11821 

0.0810 

0.1569 

0.5982 

0.5225 

0.4733 

0.6304 

0.6277 

44 

44/WI/0/9/12 

10288 

0.0929 

0.2148 

0.5372 

0.4791 

0.4185 

0.7822 

0.6358 

45 

44/WI/1 5/3/0 

14197 

0.0935 

0.1573 

0.5974 

0.6338 

0.5566 

0.8196 

0.7988 

46 

44/WI/1 5/6/0 

11353 

0.1115 

0.2164 

0.4974 

0.5538 

04548 

0.7973 

0.8108 

47 

44/WI/1 5/9/0 

11746 

0.1173 

0.2350 

0.4631 

0.5298 

0.4133 

0.7981 

0.8101 

48 

61/SU/3/3/12 

11025 

0.0843 

0.1557 

0.6185 

0.6461 

0.5790 

0.8107 

0.8639 

o 

61/SU/3/6/12 

10020 

0.1155 

0.2057 

0.5250 

0.5580 

0.5247 

0.7676 

0.8651 

50 

61/SU/3«/12 

8977 

0.1156 

0.2070 

0.5347 

0.5579 

0.5307 

0.7716 

0.6635 

El 

61/SU/l  5/3/0 

13866 

0.0949 

0.1645 

0.5692 

0.6310 

0.5107 

0.8307 

0.8095 

m 

61/SU/15/6/0 

12354 

0.1170 

0.2130 

0.4803 

0.5717 

0.4175 

0.8621 

0.8095 

M 

61/SU/l  5/9/0 

12014 

0.1165 

0.2131 

0.4797 

0.5702 

0.4148 

0.8668 

0.8114 

Q 

61/WI/3/3/12 

11043 

0.0871 

0.1381 

0.6027 

0.6653 

0.55 1 9 

0.8183 

0.8522 

55 

61/WI/3/6/12 

10166 

0.1168 

0.2111 

0.4777 

0.5461 

0.4339 

0.7561 

0.8506 

El 

61/WI/3/9/12 

9000 

0.1203 

0.2287 

0.4590 

0.5253 

0.4127 

0.7541 

0.8511 

57 

61/WI/15/3/0 

12771 

0.1062 

0.1780 

04973 

0.5917 

0.3999 

0.8235 

0.7691 

58 

61/WI/15/6/0 

11827 

0.1260 

0.2220 

0.4032 

0.5383 

0.2894 

0.8484 

0.7731 

59 

61/WI/15«/0 

10536 

0.1275 

0.2271 

0.3942 

0.5296 

0.2812 

0.8492 

0.7738 

o 

30/SU/21/3/12L 

11143 

0.0783 

0.1147 

0.6178 

0.6091 

0.5147 

0.6557 

0.7629 

61 

30/SU/21/6/0L 

10329 

0.0985 

0.1628 

0.5183 

0.4997 

0.4069 

0.5692 

0.7633 

62 

30/SU/21/9/0L 

9289 

0.1061 

0.1873 

0.4875 

0.4596 

0.3701 

0.5550 

0.7615 

1^ 

30/SU/21/3/12W 

2617 

0.0487 

0.0793 

0.7746 

0.7658 

0.7318 

0.8345 

0.8785 

El 

30/SU/21/6/0W 

2102 

0.0701 

0.1435 

0.6974 

0.6556 

0.6637 

0.7897 

0.8739 

30/SU/21/9/0W 

1719 

0.0738 

0.1573 

0.6766 

0.6358 

0.6492 

0.7865 

0.8802 

66 

61/WI/3/3/12L 

8331 

0.0883 

0.1397 

0.5942 

0.6578 

0.5337 

0.8126 

0.8548 

67 

61/WI/3/6/12L 

7801 

0.1212 

0.2023 

0.4583 

0.5533 

0.4073 

0.7586 

0.8523 

61/WI/3«/12L 

7073 

0.1241 

0.2178 

0.4383 

0.5395 

0.3863 

0.7645 

0.8517 

69 

61/WI/3/3/12W 

1909 

0.0763 

0.1371 

0.6469 

0.6867 

0.6097 

0.8481 

0.8523 

70 

61/WI/3/6/12W 

1651 

0.0887 

0.2532 

0.5942 

0.5257 

0.5584 

0.7953 

0.8474 

El 

61/WI/3/9/12W 

1738 

0.0908 

0.3192 

0.5978 

0.4465 

0.5748 

0.7658 

0.8550 

In  this  effort  we  did  not  have  an  absolute  standard  by  which  to  evaluate  the  MOS 
equations.  Instead,  a  useful  measure  of  performance  was  the  relative  improvement  of 
MOS-based  forecasts  over  persistence  using  the  Brier  score.  In  general,  we  found  that 
Brier  scores  for  MOS-based  forecasts  were  35-55%  smaller  than  their  persistence-based 
counterparts.  A  quick  glance  at  Figures  50  through  52  shows  the  qualitative  improvement 
of  the  MOS-based  forecasts  over  persistence. 
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In  Figure  54  we  show  the  averaged  percentage  improvement  of  MOS  over  persis¬ 
tence  as  a  function  of  forecast  length  for  each  box.  That  figure  shows  that  the  relative  im¬ 
provement  of  MOS  over  persistence  increases  with  increasing  forecast  length  as  the 
influence  of  persistence  in  the  forecasts  decreases. 

The  20/20  Score 

The  20/20  score  (S20/20)  measures  the  fraction  of  the  counts  in  the  contingency  table 
that  are  within  20%  (i.e.,  4  bins)  of  the  diagonal  on  either  side.  It  is  calculated  from  the 
contingency  table  from  the  following: 

21  min(21j  +  4) 

S20/20  =  H  X  Z  (3-19) 

j  =  l  i=max(lj-4) 

The  20/20  scores  for  all  of  the  71  cases  (grouped  by  RTNEPH  box  number)  are  dis¬ 
played  in  Figures  55  through  57.  The  scores  are  also  listed  by  case  number  in  Table  15.  A 
perfect  forecast  would  have  a  20/20  score  equal  to  one,  since  all  the  data  points  would  lie 
along  the  diagonal.  A  score  of  zero  (the  minimum  possible  score)  would  correspond  to  all  of 
the  data  points  lying  outside  20%  of  the  diagonal  of  the  contingency  table. 

UnUke  the  results  from  our  analysis  of  the  Brier  score,  MOS  did  not  outperform  per¬ 
sistence  consistently  with  respect  to  the  20/20  score.  Similar  to  the  Brier  score  results  how¬ 
ever,  we  found  a  trend  toward  decreasing  skill  with  increasing  forecast  length,  as  expected. 

Results  varied  considerably  as  functions  of  the  region,  season,  and  time  of  day.  For 
example,  in  the  box  44  cases,  MOS  outperformed  persistence  (i.e.,  MOS  had  higher  20/20 
scores  than  persistence)  consistently  for  the  OZ  forecast  valid  time  (both  winter  and  sum¬ 
mer),  but  did  not  perform  as  well  for  the  15Z  forecast  valid  time.  In  box  61,  we  found  that 
MOS  outperformed  persistence  only  in  two  cases,  both  of  which  contained  only  water  grid 
points  (case  numbers  70  and  71).  In  box  30,  overall  scores  (both  persistence  and  MOS) 
were  higher  for  the  cases  with  valid  forecast  time  of  9Z  than  21Z  in  all  four  seasons  and 
were  highest  for  the  all  water  case  (summer,  21Z,  case  number  63).  S20/20  for  the  MOS- 
based  forecasts  were  higher  than  those  corresponding  to  persistence  for  many  of  the  9Z 
cases  in  box  30. 

We  found  that  for  most  box/season/times,  MOS  forecasts  were  consistently  higher 
or  lower  than  persistence  over  all  forecast  lengths.  For  example,  if  MOS  performed  better 
than  persistence  for  a  3-hour  forecast  for  a  given  box,  season,  and  valid  time,  it  generally 
performed  better  for  the  6-  and  9-hour  forecasts  also. 
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As  we  saw  in  our  analysis  of  the  Brier  score  results  of  those  Box  30  cases  in  which 
we  developed  two  forecast  models  (one  using  the  OZ  GSM  data  and  the  other  using  the  12Z 
data),  there  was  no  evidence  to  suggest  that  there  was  any  measurable  difference  in  per¬ 
formance.  That  is,  the  MOS-based  forecast  models  accounted  for  any  existing  GSM  spin- 
up  problems  as  measured  by  the  20/20  skill  score  as  well  as  the  Brier  score. 

Sharpness 

Sharpness  is  a  measure  of  a  single  cloud  cover  distribution.  It  measures  the  weight 
of  the  distribution  within  20%  (i.e.,  4  bins)  of  the  outer  bins.  That  is,  it  is  the  number  of 
data  points  whose  cloud  cover  is  less  than  22.5%  or  greater  than  77.5%. 

We  computed  the  sharpness  of  the  observed  and  forecast  (MOS-based  and  j)ersis- 
tence)  marginal  distributions  directly  from  the  contingency  table  using  the  following: 


N 


21 

I 

21 

+  S 

j=l 

i=l 

i=17 

21 

5 

21 

s 

+  Z 

i=l 

j=l 

i  =  17 

(3-20) 


(3-21) 


Sharpness  ranges  from  a  maximum  of  one,  corresponding  to  a  distribution  in  which  all  of 
the  data  points  are  in  the  outer  20%  bins,  to  a  minimum  of  zero,  corresponding  to  all  the 
data  points  falling  in  the  interior  60%  bins. 

Because  sharpness  characterizes  a  cloud  cover  distribution,  we  gain  meaning  not 
from  the  value  itself,  but  from  comparing  the  sharpness  of  one  distribution  to  another.  We 
computed  the  sharpness  for  the  ground  truth  (observed)  cloud  cover,  the  persisted  cloud 
cover,  and  the  MOS-based  forecast  cloud  cover.  Those  values  are  included  in  Table  15.  The 
bar  graphs  in  Figure  58  through  60  show  the  sharpness  values  for  observed  and  forecast 
cloud  cover  for  the  71  cases  grouped  by  RTNEPH  box  number. 

Observations  of  cloud  cover  tend  to  have  U-,  J-,  or  L-shaped  distributions  with  most 
values  occurring  in  the  outer  (clear  and  overcast)  bins.  An  example  of  this  type  of  distribu¬ 
tion  was  presented  in  Figure  12.  Because  persisted  cloud  cover  is  simply  observed  cloud 
cover  that  has  been  shifted  in  time,  we  expect  similar  sharpness  values  for  those  distribu¬ 
tions.  The  MOS-based  forecasts  were  less  sharp  than  either  the  observed  or  persisted 
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cloud  cover.  That  was  to  be  expected  since  the  MOS  technique  (based  on  least  squares 
minimization  of  errors)  tends  to  predict  conservative  (i.e.,  tending  toward  the  mean)  cloud 
cover  amounts.  Figures  61  and  62  show  observed  and  forecast  cloud  cover  distributions, 
respectively,  with  very  different  sharpness  values. 

In  fact,  MOS-based  forecasts  become  increasingly  conservative  (i.e.,  less  sharp) 
with  forecast  length,  such  that  a  3-hour  forecast  is  more  sharp  than  a  6-hour  forecast 
which  is  more  sharp  than  a  9-hour  forecast.  This  effect  is  shown  clearly  in  Figure  63  and 
has  been  reported  by  others  as  well  (Ref.  17). 

One  method  of  producing  sharper  forecasts  is  to  replace  the  least  squares  mini¬ 
mization  in  the  linear  regression  with  a  method  that  minimizes  absolute  error  instead. 
Also,  as  discussed  in  Ref.  26,  one  could  use  the  logistic  response  function  in  categorical 
cloud  amount  forecasts  to  better  model  the  shape  of  the  cloud  cover  distributions.  Both  of 
these  methods  are  more  difficult  and  expensive  than  linear  regression  because  they  re¬ 
quire  iteration  and  good  initial  estimates  of  model  coefficients. 

The  technique  that  we  employed  in  this  feasibility  study  emphasized  reduction  in 
total  squared  error  and  thus  performs  very  well  with  respect  to  the  Brier  score.  In  some 
operational  settings,  forecast  sharpness  may  be  crucial,  and  therefore  a  non-linear  tech¬ 
nique  that  yields  sharper  forecasts  would  be  preferred. 

3.5.2  Estimation  of  Performance  Using  the  Jackknife  Technique 

Tb  answer  the  question  “How  will  the  models  developed  under  this  effort  perform 
in  the  futme  with  independent  weather  data?”  we  used  the  jackknife  resampling  tech¬ 
nique  to  estimate  the  performance  of  two  models  and  determine  the  confidence  limits  on 
those  estimates. 

The  jackknife  technique  is  one  of  a  number  of  resampling  techniques  that  esti¬ 
mates  model  performance  (Ref.  12).  It  improves  upon  earlier  techniques  such  as  cross  val¬ 
idation  which  divides  the  data  set  in  two  (using  one  half  for  model  development  and  the 
other  half  for  validation)  by  using  more  of  the  data  for  development  and  computing  multi¬ 
ple  performance  measures  using  different  subsamples  of  the  data.  Due  to  the  limited  size 
of  our  data  set,  it  was  necessary  to  keep  as  much  of  the  data  as  possible  during  evaluation; 
the  jackknife  technique  allowed  us  to  do  that.  In  addition,  the  jackknife  technique  ac¬ 
counts  for  any  bias  that  might  be  present  when  small  samples  are  used  in  the  estimation. 
The  method  consists  of  the  following  steps. 
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STEP  1  Split  the  model  development  data  set  into  M  groups  of  equal  size.  We 
used  9  groups  of  10  contiguous  days  each  over  a  season. 

STEP  2  Leave  out  one  group  of  data  and  develop  the  regression  model  (i.e.,  deter¬ 
mine  model  coefficients).  Repeat  the  process  for  each  of  the  M  groups, 
leaving  out  a  different  group  of  data  in  each  model. 

STEP  3  Validate  each  model  with  the  omitted  group  of  data.  Compute  a  perfor¬ 
mance  statistic.  We  used  the  Brier  score. 

STEP  4  Use  the  sample  of  M  performance  statistics  along  with  the  performance 
computed  using  the  whole  data  set  to  compute  average  performance  and 
the  variance  of  the  average. 

We  selected  two  model  development  data  sets  from  the  group  of  71  for  which  to  per¬ 
form  the  jackknife  analysis.  The  first  case  (#12  from  box  30)  had  the  lowest  Brier  score  (i.e., 
the  highest  performance)  of  all  the  71  cases.  The  second  case  (#59  frY>m  box  61)  had  the  high¬ 
est  Brier  score  (i.e.,  lowest  performance)  of  all  the  cases.  It  also  contained  fewer  grid  points 
due  to  the  geometry  of  the  box  and  thus  the  presence  of  off-world  points.  We  selected  these 
two  cases  to  provide  upper  and  lower  bounds  for  the  results  of  our  jackknife  analysis. 

The  jackknife  estimates  are  computed  from  the  M  individual  performance  esti¬ 
mates  as  follows: 

yaii  is  the  estimate  using  all  data 

yj  is  the  estimate  with  the  jth  group  omitted 

yj*  is  a  “pseudo”  estimate  defined  by 

yj  =  M  yai,  -  (M-1)  y^  (3-22) 

y*  is  the  jackknife  performance  estimate 

M 

y*  =  M  ^3-23) 

J=i 

s?  is  the  estimated  variance  of  y* 

""  mu  - 1)  Z  fyj  “  y  0-24) 

j=i 
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The  yaii,  yj,  y*,  and  s*  values,  computed  for  the  two  cases,  are  presented  in  Table  16. 
There  is  little  change  in  estimated  performance  using  the  jackknife  analysis  (y*)  versus 
using  the  whole  data  set  (yaii)  in  both  cases.  This  implies  that  the  MOS  models  are  robust 
with  respect  to  the  variability  in  the  data  sets. 

In  a  preliminary  analysis  we  divided  the  data  set  into  9  groups  by  sampling  every 
ninth  grid  point  from  the  15  x  15  grid  points  in  the  RTNEPH  box  over  one  season.  (Recall 
that  we  used  15  x  15  grid  points  to  allow  us  to  compute  horizontal  gradients  for  diver¬ 
gence,  vorticity,  etc.)  We  found  that  the  variance  in  the  estimated  performance  was  about 
a  factor  of  five  smaller  than  the  numbers  shown  in  Table  16  which  were  generated  using 
the  9  10-day  data  sets.  This  was  due  to  the  fact  that  our  preliminary  sampling  technique 
which  included  every  ninth  g^d  point  had  not  accounted  for  the  high  spatial  correlation 
in  the  cloud  cover  data.  Thus  our  9  model  data  sets  were  highly  correlated.  By  sampling 
9  10-day  periods  instead,  as  we  did  in  our  final  analysis,  the  model  development  data  sets 
were  less  correlated.  However,  due  to  the  nature  of  cloud  cover  data,  the  resulting  vari¬ 
ances  (those  reported  in  Table  16)  are  most  likely  still  overly  optimistic.  Ideally,  one  would 
evaluate  the  performance  of  the  MOS-based  models  against  completely  independent  data 
taken  from  another  year. 


Table  16  Results  of  Jackknife  Analysis  For  Cases  #12  and  #59.  The  Average 
Performance  (y*)  and  Variance  ( sf )  for  Both  Cases  are  Highlighted 


Case  #12 

Case  #59 

Vail 

.0589 

.1275 

yi 

.0614 

.1342 

y2 

.0595 

.1308 

ya 

.0643 

.1326 

y4 

.0658 

.1134 

ys 

.0551 

.1266 

ye 

.0577 

.1419 

y? 

.0544 

.1369 

ye 

.0549 

.1313 

yg 

.0586 

.1263 

• 

y 

.0575 

.1039 

si 

5.3320-4 

4.5780-4 
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3.5^  Verifying  MOS-based  Probabilities  of  Cloud  Cover  Cat^ories 

We  used  a  Bayesian  approach  to  assess  the  likelihood  that  the  MOS-based  forecast 
probabilities  are  “correct”.  That  is,  that  they  accurately  represent  the  frequency  of  occur¬ 
rence,  in  nature,  of  the  events  of  interest.  For  this  evaluation,  we  selected  the  same  two 
model  development  data  sets  used  in  the  evaluation  of  the  MOS-based  total  cloud  amount 
forecasts  presented  in  Section  3.5.:  Case  #12  from  box  30  and  Case  #59  from  box  61.  These 
two  cases  represent  the  highest  and  lowest  performance  (as  measured  by  the  accuracy  of 
the  total  cloud  amount  forecasts)  out  of  the  set  of  71  listed  in  Tbble  3. 

The  following  paragraphs  describe  our  approach  in  general  terms  and  the  way  in 
which  the  available  data  are  used  in  its  implementation.  Following  that  description,  a  key 
formula  based  on  the  Bayesian  viewpoint  is  presented.  (Derivation  of  this  formula  is  pro¬ 
vided  in  Appendix  B.)  Finally,  the  results  are  presented  in  tabular  and  graphical  forms 
accompanied  by  a  discussion  of  their  implications. 

The  Bayesian  Viewpoint 

Tb  carry  out  the  evaluation  task,  we  must  provide  a  careful  definition  of  the  “events 
of  interest”  mentioned  above.  MOS-based  probability  forecasts  have  been  developed  for 
71  cases,  each  defined  by  a  region  (i.e.,  an  RTNEPH  box),  a  season,  a  forecast  valid  time, 
and  a  forecast  length.  In  each  of  the  two  cases  we  analyzed,  we  considered  27  “events  of 
interest”  as  follows: 

For  each  one  of  nine  cloud-cover  fraction  (ccf)  bins,  Bj,  (i.e.,  Bi  =  “0.1”  represents 

{5%  <  ccf  <=  15%}  ,  B2  =  “0.2”  represents  {15%  <  ccf  <=  25%|,  etc.) 

and 

For  each  one  of  three  specified  cloud-cover  threshold  values,  Tj,  {0.3,  0.5,  0.7}, 

The  event  E(i,  j)  occurs  at  the  specified  valid  time  if  the  MOS  cloud  cover  forecast 

for  that  time  is  in  bin  B,  and  the  actual  ccf  at  that  time  is  less  than  or  equal  to  Tj. 

For  an  ideal  validation,  one  would  observe  several  thousand  occurrences  of  each  ccf 
forecast  bin  for  the  region,  season,  valid  time  and  forecast  length  of  interest.  One  could 
then  count  the  number  of  times  that  the  actual  ccf  is  less  than  each  threshold  of  interest 
and  accept  or  reject  h3rpotheses  regarding  the  validity  of  the  MOS  forecast  probabilities 
based  on  the  observed  frequency  of  occurrence. 
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Since  such  vast  amounts  of  data  are  not  available  for  any  one  model  case,  we 
adopted  a  Bayesian  viewpoint  in  which  MOS-based  probabilities  of  the  27  events,  E(i,  j), 
are  generated  using  one-half  of  the  available  data  (Ref.  18).  These  probabilities  are  then 
viewed  as  prior  information,  establishing  a  prior  probability  measure  on  the  unknown, 
true  probability  that  E(i,  j)  will  occur.  The  other  half  of  the  original  data  is  then  used  to 
update  that  prior  probability  measure,  obtaining  a  posterior  probability  measure  on  the 
true  probability.  If  the  new  data  substantially  increase  the  weight  at  the  MOS-based  val¬ 
ue,  we  will  have  shown  consistency  among  at  least  the  limited  data  set  available.  Con¬ 
versely,  if  the  weight  is  drastically  decreased  at  the  MOS  value,  then  applicability  of  the 
MOS  model  for  the  given  forecast  valid  and  lead  times)  is  called  into  question. 

A  Formula  for  Posterior  Probability 

Results  presented  in  this  section  are  calculated  using  the  following  formula  for  the 
posterior  probability  (p*)  that  the  true  probability  of  the  event  of  interest  (p),  is  equal  to 
the  MOS  forecast  probability  (p). 

p*  = _ (3-25) 

Up/Dg)  +  I(n,k) 

A  k  A  A 

p  (l-p)"’**  =  likelihood  of  p  given  data  set  Dg 

n  =  total  number  of  observations  in  data  set  Dg 
k  =  number  of  times  that  E(i,j)  occurred  in  D9 

and 

1 

I(n,k)  =  j  x*Hl-x)"  *^  dx.  (3-26) 

0 

To  evaluate  the  integral,  we  used  the  following  M-term  approximation: 

M 

I(n,k)  =  ^  x|^(l-x,n)"-‘^Am  (3-27) 

m=  1 


where 

Up/Dg)  = 
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where  (An,,  Xm)  are  the  width  and  midpoint  of  the  interval  in  the  sum.  All  of  the  results 
presented  here  were  calculated  with  M  =  100  and  A^  =  0.01.  IVpically,  only  fifteen  terms 
of  the  sum  surrounding  the  maximum  of  the  integrand  were  required  to  evaluate  the  sum 
accurately  to  six  significant  digits.  This  is  because  of  the  integrand’s  extremely  narrow 
peak  caused  by  the  large  values  of  n  and  k  for  most  cases  (e.g.,  n  =  729  and  k  =  568  for  case 
#12,  ccf  category  0.3,  with  threshold  =  0.5). 

Results 

We  considered  two  cases,  each  defined  by  a  region  (i.e.,  an  RTNEPH  box),  season, 
forecast  valid  time,  and  forecast  length: 

Case  #12  —  Box  30,  summer,  three  hour  forecast  valid  at  09Z 
Case  #59  —  Box  61,  winter,  nine  hour  forecast  valid  at  15Z. 

The  data  sets  for  each  case  were  divided  into  the  first  and  second  halves  of  each  season. 
Data  set  1  (columns  2, 5  and  8  of  Tables  17  and  18)  was  used  to  develop  MOS-based  models 
of  probability  as  described  in  Section  3.4.2  of  this  report.  The  “Data  Set  1”  columns  in  the 
tables  are  the  resulting  MOS-forecast  probabilities.  Data  set  2  was  then  used  as  “new” 
data  for  verifying  the  correctness  of  those  forecasts.  The  “Data  Set  2”  columns  in  the  tables 
are  the  frequencies  of  occurrence  as  described  above. 

If  the  MOS  models  are  any  good  at  all,  we  would  expect  that  the  MOS  forecast  prob¬ 
abilities  and  the  fi*equencies  of  occurrence  would  be  reasonably  close  to  each  other.  In  fact, 
looking  at  Tables  17  and  18,  this  is  generally  the  case.  But,  are  the  differences  significant? 
A  Bayesian  viewpoint  will  help  to  answer  that  question. 

We  will  discuss  two  examples,  forecast  ccf  categories  0.6  and  0.7  with  threshold  = 
0.3  in  each  example,  and  use  the  data  for  case  #59  (Table  18)  to  illustrate  the  approach. 
The  events  of  interest  are  that  the  actual  ccf  is  less  than  30%  when  the  MOS  forecast  ccf 
is  between  55%  and  65%  or  between  65%  and  75%,  respectively.  We  assume  that  the  un¬ 
known,  true  probability  is  a  random  variable  distributed  on  [0, 1]  with  50%  of  the  proba¬ 
bility  at  the  MOS- value  (0.30  and  0.18  for  our  two  examples)  and  the  other  50% 
distributed  uniformly  on  [0,  1].  We  then  calculate  the  posterior  probability  at  the  MOS- 
value  based  on  processing  the  new  (Data  Set  2)  data.  These  posterior  probabilities  are, 
from  the  columns  labeled  “Probability”  in  Table  18,  0.31  and  0.96  for  our  two  examples. 
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Table  17  MOS  Probability  Verification  (Case  #12) 
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Thus,  for  the  first  example  (MOS  forecast  category  0,6  and  a  30%  ccf  threshold)  con¬ 
sidering  the  second  data  set  reduced  our  confidence  that  the  MOS  forecast  probability  (0.30) 
is  correct  firem  0,5  (the  prior  probability)  down  to  0.31.  On  the  other  hand,  for  a  MOS  fore¬ 
cast  ccf  in  the  0.7  bin,  the  probability  is  increased  substantially  fiinm  0.5  up  to  0.96.  In  the 
former  case,  considerable  variability  in  the  data  over  the  season  of  interest  is  indicated, 
casting  doubt  on  the  “correctness”  of  the  resulting  MOS  model.  In  the  latter  case,  however, 
the  second  data  set  considerably  strengthened  our  confidence  in  the  MOS  forecast,  indicat¬ 
ing  that  at  least  for  the  one  season  of  data  available,  the  MOS  forecast  probability  is  correct. 

In  the  two  tables  we  have  used  shading  to  identify  cases  in  which  the  MOS  forecast 
probability  is  “strongly  rejected”  and  “marginally  accepted”  (or  “neutral”).  Somewhat  ar¬ 
bitrarily,  we  used  posterior  probability  ranges  [0.0, 0.3]  and  [0.3, 0.8]  to  define  these  two 
designations. 

Figures  64  through  69  present  cross-sections  of  the  data  from  the  two  tables  which 
highlight  the  sensitivity  of  our  results  to  variations  between  the  two  parts  of  the  parti¬ 
tioned  data  sets.  In  each  figure,  the  line  graphs  show  the  MOS  forecast  probabilities  (solid 
squares,  Data  Set  1)  and  frequencies  of  occurrence  (hollow  squares.  Data  Set  2)  as  func¬ 
tions  of  the  forecast  ccf  category.  The  superimposed  bars  represent  the  posterior  probabili¬ 
ties  for  each  case  and  clearly  show  those  cases  for  which  the  data  set  partitions  are  most 
inconsistent. 

While  it  is  difficult  to  generalize  regarding  why  certain  of  the  MOS  forecasting 
models  are  rejected,  some  understanding  of  the  causes  can  be  obtained  by  examining  a  few 
specific  examples.  We  will  consider  the  case  #12  (box  30)  results  found  in  Table  17  and 
plotted  in  Figures  64  through  66  and  examine  the  evaluation  results  for  three  specific  ccf 
forecast  categories  in  more  detail. 

For  ccf  forecast  category  0.7,  Table  17  reveals  that  the  MOS  model  was  strongly  af¬ 
firmed  for  all  three  thresholds  (0.3, 0.5,  ana  o  7).  Figure  70  contains  plots  of  the  experimen¬ 
tal  distribution  functions  of  the  true  ccf  (i.e.,  the  RTNEPH  anal3rsis  values)  for  only  those 
occasions  at  which  the  MOS-based  ccf  forecast  category  was  0.7.  In  the  figure,  arrows  indi¬ 
cate  the  threshold  locations;  the  intercepts  of  the  Data  Set  1 2uid  2  curves  with  those  thresh¬ 
olds  define  the  probabilities  foimd  in  Table  17  under  Forecast  Cloud-Cover  Fraction  0.7.  For 
this  forecast  category,  the  two  distribution  functions  are  close  to  another  at  the  threshold 
values,  consistent  with  the  strong  affirmation  of  the  MOS-based  probability. 
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Contrast  the  above  case  with  that  for  ccf  forecast  category  =  0.2  whose  experimen¬ 
tal  distribution  functions  are  plotted  in  Figure  71.  In  this  case,  because  of  the  substantial¬ 
ly  larger  percentage  of  occxirrences  of  clear  conditions  (ccf  <  0.15)  in  Data  Set  2  (66%)  than 
Data  Set  1  (60%),  the  distribution  intercepts  with  the  0.3  threshold  line  are  farther  apart 
(0.752  and  0.794  are  the  actual  vgdues).  Thus,  for  the  threshold  =  0.3  case,  the  MOS  model 
was  strongly  rejected  by  the  evaluation  procedure,  while  for  the  two  higher  threshold,  no 
inconsistency  was  noted. 

Similarly,  for  ccf  forecast  category  =  0.8,  shown  in  Figure  72,  the  distributions  are 
rather  far  apart  at  the  0.7  threshold  (and  the  MOS  model  is  strongly  rejected)  while  the 
models  are  supported  for  the  two  lower  threshold  values.  In  this  case,  the  discrepancies 
arise  because  Data  Set  2  contains  a  significantly  larger  percentage  of  cloudy  conditions 
(69%)  than  Data  Set  1  (59%). 

We  can  summarize  the  above  as  follows:  the  MOS  probability  forecasts  were  re¬ 
jected  when  the  differences  between  the  cloud-cover  fraction  distributions  of  the  original 
and  verifying  data  sets  occurred  at  ccf  values  near  the  threshold  value  of  interest.  (I.e., 
differences  in  the  fraction  of  clear  cases  caused  the  threshold  =  0.3  model  to  be  rejected, 
and  differences  in  the  fraction  of  very  cloudy  cases  cause  the  threshold  =  0.7  model  to  be 
rejected.)  Basing  the  MOS  model  on  the  complete  data  set  (instead  of  on  the  earliest  50% 
as  in  this  evaluation)  should  improve  the  model  robustness  with  respect  to  seasonal  varia¬ 
tions.  Whether  or  not  that  improvement  will  be  sufficient  must  be  evaluated  using  the 
multi-yeeir  data  sets. 


OVERALL  CONCLUSIONS 

1.  For  most  events  (39  out  of  the  54  examined),  the  MOS  forecast  values  are  strongly 
affirmed  by  the  second  data  set.  This  indicates,  at  least,  that  within  the  one  season 
for  which  we  have  data,  the  MOS-models  for  these  two  cases  are  reasonably  good. 

2.  For  5  out  of  the  54  cases,  the  MOS  forecast  values  were  strongly  rejected.  However, 
comparing  the  line  plots  in  Figures  64  through  69  shows  that,  even  for  those  five 
cases,  the  differences  between  the  fi:«quen(y  of  occurrence  and  the  MOS  forecast  are 
not  very  large.  Some  variation  within  the  season  is  indicated,  but  the  MOS  models 
are  still  not  drastically  in  error.  (For  example;  for  case  #12,  box  30,  with  MOS  forecast 
ccf  =  0.2,  the  forecast  probability  that  the  actual  ccf  would  be  less  than  0.3  was  0.75, 
but  in  fact,  in  Data  Set  2,  79%  of  the  cases  had  actual  ccf  less  than  0.3.) 
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Figure  26  Example  ofan  Error  Histogram  Showing  How  to 

Compute  Conditional  Cloud  Cover  Probability 


BOX  38,  WINTER 
Eigenvalues 


PCI  PC2  PC3  PC4  PCS 


Figure  27  Variance  of  Principal  Components  (PC)  Computed 

For  Box  38,  Winter 
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BOX  38,  WINTER 
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Figure  29  Breakdown  of  the  a)  First  and  b)  Second  Principal  Component 

Into  its  Component  Parts,  Box  38,  Winter 
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BOX  61,  WINTER  BOX  61,  WINTER 
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Figure  30  Breakdown  of  the  a)  First,  b)  Second,  and  c)  Third  Principal 
Component  Into  its  Component  Parts,  Box  61,  Winter 
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FtforeSS  Schematic  Showing  the  Grid  Geometiy  Used  to  Compute 

Derived  Variables  Such  as  Divergence  and  Vorticity.  The 
(x,y)  Coordinate  System  Used  in  this  Figure  Corresponds  to 
the  (i  J)  Coordinate  System  Used  in  the  Nephanalysis  Fields 
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FigoreSd  Definition  ofLow,  Middle,  and  Hi^  Layers  in  Terms  of 

Sigma  Surfaces  Used  in  Model  Development  (Ref  24) 
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Figure  35  Histogram  Showing  All  Those  Predictors  That  Were  Selected  in 

the  Top  15  More  Than  25%  of  the  Time  Over  all  7 1  Models 
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Normal  Probability  Plot  of  Residuals 
Box  30,  Fall,  21Z,  3-Hour  Forecast 


Residuals 

Figure  36  Normal  Plot  ofthe  Residuals  for  Case  #24 
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Expected  Normal  Value 


Normal  Probability  Plot  of  Residuals 
Box  30,  Fall,  21Z,  9-Hour  Forecast 
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Figure  38  Normal  Plot  ofthe  Residuals  for  Case  #27 
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Expected  Normal  Value 


Normal  ProbcUoility  Plot  of  Residuals 
Box  61,  Summer,  15Z,  6-Hour  Forecast 


Residuals 

Figure  40  Normal  Plot  ofthe  Residuals  for  Case  #52 

Showing  a  Poor  Fit  to  the  “Normal  Line” 
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Figure  44  MOS-based  Probability  Forecast  That  the  Observed  Cloud  Cover  Will 

be  Less  Than  or  Equal  to  the  Specified  Cloud  Cover  Threshold  Given 
the  MOS-predicted  Cloud  Amount  (cc)  in  Case  #24 
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Figure  45  MOS-based  Probability  Forecast  That  the  Observed  Cloud  Cover  Will 

be  Less  Than  or  Equal  to  the  Specified  Cloud  Cover  Threshold  Given 
the  MOS-predicted  Cloud  Amount  (cc)  in  Case  #25 
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Figure  46  MOS-based  Probability  Forecast  That  the  Observed  Cloud  Cover  Will 

be  Less  Than  or  Equal  to  the  Specified  Cloud  Cover  Threshold  Given 
the  MOS-predicted  Cloud  Amount  (cc)  in  Case  #27 


114 


Cloud  Cover  Threshold  (%) 


Figure  47  MOS-based  Probability  Forecast  That  the  Observed  Cloud  Cover  Will 

be  Less  Than  or  Equal  to  the  Specified  Cloud  Cover  Threshold  Given 
the  MOS-predicted  Cloud  Amount  (cc)  in  Case  #42 
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Figure  48  MOS-based  Probability  Forecast  That  the  Observed  Cloud 

be  Less  Than  or  Equal  to  the  Specified  Cloud  Cover  Threshold 
Given  the  MOS-predicted  Cloud  Amount  (cc)  in  Case  #52 
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Forecast  Cloud  Cover 
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Figure  5 1  Comparison  of  Brier  Scores  for  MOS-based  and  Persistence-based  Forecasts  for  All  Cases  in  Box  44 
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Case  Number 

Figure  52  Comparison  of  Brier  Scores  for  MOS-based  and  Persistence-based  Forecasts  for  All  Cases  in  Box  6 1 
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Figure  54  Averaged  Percentage  Improvement  in  Brier  Scores  of  MOS 

Forecasts  Over  Persistence  by  Box  and  Forecast  Length 
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Figure  65  Comparison  of 20/20  Scores  for  MOS-based  and  Persistence-based  Forecasts  for  all  Cases  in  Box  30 
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Figure  56  Comparison  of 20/20  Scores  for  MOS-based  and  Persistence-based  Forecasts  for  all  Cases  in  Box  44 


Figure  57  Comparison  of 20/20  Scores  for  MOS-based  and  Persistence-based  Forecasts  for  all  Cases  in  Box  6 1 
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Figure  58  Comparison  of  Sharpness  for  MOS-based  Forecasts,  Persistence-based  Forecasts  and  Observed 

Cloud  Cover  for  all  Cases  in  Box  30 
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Figure  69  Comparison  of  Sharpness  for  MOS-based  Forecasts,  Persistence-based  Forecasts  and  Observed 

Cloud  Cover  for  all  Cases  in  Box  44 
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Figure  60  Comparison  of  Sharpness  for  MOS-based  Forecasts,  Persistence-based  Forecasts  and  Observed 

Cloud  Cover  for  all  Cases  in  Box  6 1 
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Figure  61  Histogram  of  Observed  Cloud  Cover  Values  for  Case  #32  (Box  30, 
Winter,  2 IZ)  Showing  a  Very  Sharp  U-shaped  Distribution 
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MOS  Forecast  Cloud  Cover  (%) 

Figure  62  Histogram  of  3>Hour  Forecast  Cloud  Cover  Values  for  Case  #32  (Box 

30,  Winter,  21Z)  Showing  a  Much  Smoother  Distribution  Than  the 
Observed  Cloud  Cover  Shown  in  Figure  61 
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Figure  63  Averaged  Sharpness  for  Each  of  the  Box  Regions  Studied  as 

a  Function  of  Forecast  Length 
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Figure  64  Verifying  Forecast  Cloud  Cover  Probability  Case  #  12  Cloud  Cover  Threshold  30% 
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Figure  65  Verifying  Forecast  Cloud  Cover  Probability  Case  #12  Cloud  Cover  Threshold  60% 
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Figure  66  Verifying  Forecast  Cloud  Cover  Pn^ability  Case  #  12  Cloud  Cover  Threshold  70% 
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Figure  67  Verifying  Forecast  Cloud  Cover  Probability  Case  #59  Cloud  Cover  Threshold  30% 


Figure  68  Verifying  Forecast  Cloud  Cover  Probability  Case  #59  Cloud  Cover  Threshold  50% 
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Figure  69  Verifying  Forecast  Cloud  Cover  Probability  Case  #59  Cloud  Cover  Threshold  70% 


Pr{  CCF  <=  ccf  }  —  Forecast  CCF  Category:  0.7 


AiiiiqeqoJd 


138 


Figure  70  Experimental  Distribution  Functions  of  the  TVue  ccf  for  the  MOS  ccf  Forecast  Category  of  0.7 
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Pr{  CCF  <=  ccf  }  —  Forecast  CCF  Category:  0.8 
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4. 


RECOMMENDATIONS 


As  discussed  in  the  previous  section,  the  MOS-based  forecast  models  developed 
under  this  effort  performed  well  in  terms  of  overall  variance  reduction  and  robustness. 
GSM  variables  added  substantially  to  the  skill  of  the  cloud  cover  forecasts  (we  saw  35-55% 
improvement  when  compared  with  persistence  forecasts  alone).  Our  analysis  showed  that 
the  MOS  approach  using  Global  Spectral  Model  variables  in  addition  to  persistence  shows 
promise  as  a  low  cost  cloud  forecasting  tool.  (We  use  the  term  "low  cost”  because  the  MOS 
approach  uses  existing  and  previously  verified  forecast  models  along  with  existing  cloud 
analysis  databases.  In  addition,  the  MOS  model  as  described  here  uses  linear  regression 
techniques  that  are  computationally  efficient.) 

However,  throughout  this  report  we  have  also  pointed  out  the  limits  of  the  models 
developed  under  this  effort.  In  this  section,  we  discuss  those  limitations  and  offer  recom¬ 
mendations  in  three  general  areas:  model  development  data,  the  modeling  approach,  and 
performance  analysis  and  validation. 


4.1  MODEL  DEVELOPMENT  DATA 

Develop  Models  Using  Higher-Resolution  GSM  Data 

We  generated  MOS-based  forecasts  of  total  cloud  cover  at  eighth-mesh  grid  resolu¬ 
tion.  Potential  predictors  at  that  resolution  included  terrain  type,  elevation,  and  persisted 
cloud  amount.  The  remaining  predictors  used  during  development  were  derived  fix>m  the 
40-wave  Global  Spectral  Model  with  a  grid  resolution  over  four  times  more  coarse  than 
eighth-mesh.  Cloud  processes  important  at  eighth-mesh  resolution  (such  as  local  terrain  ef¬ 
fects,  sea  breeze  effects,  and  localized  convection)  are  not  modeled  with  the  40-wave  (xSM. 
Therefore  the  MOS  equations  developed  with  the  40-wave  model  data  cannot  support  cloud 
formation  and  evolution  at  that  resolution.  In  addition,  the  2.5°  grid  used  in  the  40-wave 
spectral  model  effectively  contains  area  averages  of  the  model  variables  with  respect  to  the 
much  finer-resolution  eighth-mesh  grid,  resulting  in  forecasts  that  are  too  smooth. 

We  recommend  using  higher-resolution  GSM  data  in  model  development  such  as 
is  currently  available  with  the  80-wave  (and  higher)  versions  of  the  spectral  model.  The 
recently  released  Eta  model  from  the  National  Weather  Service  is  a  good  example  of  the 
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recent  advances  in  global  weather  modeling;  it  has  a  horizontal  grid  resolution  of  approxi¬ 
mately  80  km  with  38  levels  in  the  vertical.  In  the  future,  even  higher-fidelity  models  will 
be  available.  Of  course,  higher  model  resolution  alone  is  not  sufficient  to  guarantee  im¬ 
proved  MOS  forecasts.  The  models  themselves  must  accurately  predict  the  growth  and 
decay  of  disturbances  at  those  resolutions. 

The  increased  spatial  resolution  would  not  necessarily  increase  data  processing 
and  storage  costs  if  a  similar  number  of  grid  point  locations  (such  as  the  15  x  15  grid  points 
per  RTNEPH  box  used  in  this  analysis)  were  modeled.  However,  increases  in  the  number 
of  vertical  levels  can  quickly  increase  the  number  of  potential  predictors.  If  necessary,  the 
number  of  variables  in  the  vertical  could  be  reduced  by  using  only  their  most  significant 
principal  components. 

As  with  any  significant  change  to  the  driving  model,  at  least  two  years  of  stable 
development  data  would  be  required  to  develop  operational  MOS  equations  based  on 
higher-resolution  model  output. 

Include  RTNEPH  Layer  Cloud  Amounts  and  Types 

Cloud  layer  fields  reported  in  the  RTNEPH  analysis  grids  (or  available  from  pilot 
reports,  upper  air  analyses,  etc.)  could  be  used  in  conjunction  with  the  humidity  profiles 
generated  in  the  GSM.  Using  cloud  layer  data  would  allow  us  to  adjust  humidity  layers  for 
the  thickness  of  the  cloud  layer  within  each  sigma  level  and  thus  eliminate  any  ambiguity 
concerning  whether  or  not  a  cloud  layer  straddles  two  or  more  sigma  levels.  This  would  re¬ 
sult  in  better  estimi  es  for  the  maximum  layer  humidity,  a  strong  predictor  for  cloud 
amount.  The  quality  of  the  cloud  layer  data,  the  additional  value  it  might  bring  to  the  MOS 
equations,  and  the  increased  data  processing  and  storage  costs  would  all  need  to  be  weighed 
before  deciding  whether  or  not  to  include  these  data  in  future  modeling  efforts. 

Use  Independent  Persistence  and  Response  Data 

In  this  effort  we  used  the  RTNEPH  analysis  grids  both  as  persisted  cloud  cover  and 
“ground  truth.”  Thus,  there  is  potentiaUy  a  high  correlation  between  the  two.  We  at¬ 
tempted  to  remove  the  temporal  correlation  by  removing  all  data  points  which  had  not 
been  updated  in  the  three  hours  prior  to  the  forecast  valid  time.  However,  the  time  flag 
available  to  verify  data  age  only  indicates  the  age  of  the  newest  data  source  for  the  grid 
point;  not  necessarily  the  age  of  the  total  cloud  cover  reported  for  that  grid  point.  There¬ 
fore,  some  old  data  may  still  be  present  in  the  development  data  sets. 
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In  considering  the  results  of  this  feasibility  study,  it  is  important  to  note  the  relative 
improvement  of  the  MOS  approach  over  and  above  persistence  to  account  for  potential 
artificial  correlations  in  the  development  data  sets  due  to  the  age  of  the  data.  In  order  to 
assess  the  absolute  impact  of  persistence  and  MOS*based  forecasts,  independent  persis¬ 
tence  and  response  data  sources  should  be  used.  Candidate  data  sources  other  than  the 
RTNEPH  analysis  fields  include:  surface  observations  and  GOES-based  cloud  analyses. 

4.2  MODELING  APPROACH 

Include  Additional  Years  of  Data 

It  is  well  known  that  meteorological  conditions  vary  significantly  from  year  to  year. 
We  used  one  year  of  data  to  develop  MOS  models  in  this  feasibility  study.  Although  accept¬ 
able  for  a  feasibility  study  such  as  this,  it  would  be  imperative  to  use  additional  years  of 
data  (at  least  two)  for  operational  use.  Once  the  modeling  mechanism  is  defined,  process¬ 
ing  additional  data  is  a  mechanical  process.  By  using  additional  years  of  data  we  could 
achieve  more  stable  regression  coefficients  along  with  a  predictor  set  that  is  more  mean¬ 
ingful  and  less  subject  to  any  single  year’s  meteorological  peculiarities. 

Minimize  Absolute  Error 

Linear  multiple  regression  algorithms,  like  the  one  used  in  this  feasibility  study, 
use  least  squares  techniques  to  determine  the  MOS  equations.  This  accentuates  those 
cases  with  large  residuals  and  can  result  in  an  overly  conservative  fit  to  the  data  (i.e.,  ov¬ 
erly  smooth  forecasts).  To  make  sharper  forecasts,  it  may  be  desirable  to  minimize  abso¬ 
lute  error  instead  of  squared  error.  The  additional  value  of  using  absolute  error  would 
need  to  be  evaluated  against  the  possibly  significant  increase  in  computational  time  to 
perform  the  regression. 

Use  Logistic  Regression 

Logistic  regression  (Ref.  1)  is  a  computationally  intensive  regression  technique.  In 
Section  3.1.1  we  noted  its  advantages  for  modeling  binary  response  functions  such  as 
categorical  cloud  amounts.  The  logistic  function  does  not  produce  responses  that  lie  out¬ 
side  the  range  [0,  1]  as  linear  regression  techniques  do.  Future  efforts  may  evaluate  the 
need  for  a  non-linear  regression  technique  such  as  logistic  regression  by  weighing  the  sig¬ 
nificant  additional  computational  costs  (logistic  regression  is  an  iterative  method  which 
requires  good  initial  estimates  of  the  model  parameters)  against  its  advantages. 
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Use  Neural  Networks 


Another  non-linear  modeling  technique  that  is  being  used  more  and  more,  is  neural 
networks.  Neural  networks  provide  a  modeling  capability  that  is  different  from  conven¬ 
tional  methods.  Instead  of  using  a  particular  functional  dependence  between  the  mea¬ 
sured  data  and  desired  prediction,  neural  networks  adjust  themselves  to  give  an  optimal 
(with  respect  to  some  metric,  such  as  the  mean  square  error)  non-linear  transformation 
between  predictors  and  the  response  variable.  In  addition,  neural  networks  can  be  used 
to  determine  the  minimum  set  of  variables  necessary  to  solve  a  problem. 

The  availability  of  COTS  neural  net  software  tools  provides  for  a  low  risk/potentially 
high-payoff  approach  to  the  area  of  cloud  cover  prediction.  One  such  software  tool 
(ExploreNet  from  HNC)  has  already  shown  promise  for  this  problem  in  a  small  test  case 
performed  at  TASC. 


4.3  PERFORMANCE  ANALYSIS  AND  VALIDATION 

Compute  Other  Performance  Measures 

In  this  effort  we  selected  three  performance  statistics  (Sharpness,  Brier  score,  and 
20/20  score)  to  measure  the  skill  of  the  MOS-based  forecasts.  These  three  statistics  were 
selected  in  consultation  with  PL  personnel  as  meaningful  in  the  context  of  potential  ap¬ 
plications  for  this  MOS  approach.  Other  applications  may  emphasize  different  aspects  of 
the  model  forecasts.  For  example,  certain  satellite  systems  may  place  more  emphasis  on 
reliability  measures.  Future  analyses  could  include  more  or  different  performance  mea¬ 
sures  to  accommodate  the  interests  of  the  customer  community. 

Validate  More  Cases 

We  selected  two  cases  for  evaluation  from  the  71  development  cases  analyzed  in 
this  study.  Those  two  cases  were  chosen  because  they  had  the  highest  and  lowest  perfor¬ 
mance  of  all  71  cases  with  respect  to  overall  reduction  in  variance.  It  was  our  hope  that 
by  studying  these  two  extreme  cases  we  could  bound  our  error  estimates  for  both  the  total 
and  categorical  cloud  amount  forecasts.  Additional  cases  should  be  analyzed  in  the  future 
to  get  a  better  handle  on  the  variation  of  the  forecast  strengths  and  weaknesses  as  a  func¬ 
tion  of  time  of  day,  season,  and  region. 
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Validate  Against  More  Data 

In  the  two  cases  mentioned  above  we  evaluated  both  total  and  categorical  cloud 
amount  forecasts  by  developing  sample  models  using  a  portion  of  the  data  set  and  validat¬ 
ing  those  against  the  remainder  of  the  season  of  data. 

It  would  be  crucial  if  this  approach  was  to  be  taken  toward  an  operational  model, 
that  the  models  be  tested  against  additional  seasons  of  cloud  data  from  independent  years 
containing  a  greater  variety  of  weather  conditions.  The  estimates  of  model  robustness  de¬ 
termined  in  our  analysis  were  most  likely  overly  optimistic  since  the  variability  in  one 
season  is  most  likely  less  than  across  multiple  seasons.  We  recommend  developing  models 
with  at  least  two  years  of  data  and  validating  the  models  against  at  least  two  additional 
independent  years. 


APPENDIX  A 


A  list  of  all  of  the  variables  (listed  by  their  abbreviated  names)  used  in  the  pool  of 
potential  predictors  follows.  A  short  description  of  each  variable  is  provided  along  with  the 
appropriate  units  and  any  specific  numerical  ranges  if  applicable. 


VARIABLE 

NAME 

VARIABLE  DESCRIPTION 

UNITS/ RANGE 

DAY 

day  of  the  month  (1-31) 

[1.311 

TYPE 

terrain  type  at  eighth-mesh  grid  point  at  center  (lower  right)  of  half-mesh  cell 

0«  water, 

1  -  land, 
3-coast 

ELEV 

terrain  elevation  at  eighth-mesh  grid  point  at  center  (lower  right)  of  half¬ 
mesh  cell 

meters 

LAT 

latitude  at  center  of  grid  point 

degrees 

LON 

longitude  at  center  of  grid  point 

degrees 

ZENITH 

sine  of  solar  zenith 

UL 

low  level  wind  in  zonal  direction 

meters/sec 

UM 

middle  level  wind  in  zonal  direction 

meters/sec 

D3UL 

3-hr  change  in  low  level  zonal  wind 

meters/sec 

D3UM 

3-hr  change  in  mid  level  zonal  wind 

meters/sec 

06UL 

6-hr  change  in  low  level  zonal  wind 

meters/sec 

D6UM 

6-hr  change  in  mid  level  zonal  wind 

meters/sec 

09UL 

9-hr  change  in  low  level  zonal  wind 

meters/sec 

D9UM 

9-hr  change  in  mid  level  zonal  wind 

meters/sec 

VL 

low  level  wind  in  meridional  direction 

meters/sec 

VM 

middle  level  wind  in  meridional  direction 

meters/sec 

D3VL 

3-hr  change  in  low  level  meridionai  wind 

meters/sec 

D3VM 

3-hr  change  in  mid  level  meridional  wind 

meters/sec 

D6VL 

6-hr  change  in  low  level  meridional  wind 

meters/sec 

D6VM 

6-hr  change  in  mid  level  meridional  wind 

meters/sec 

D9VL 

9-hr  change  in  low  level  meridional  wind 

meters/sec 

D9VM 

9-hr  change  in  mid  level  meridional  wind 

meters/sec 
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VARIABLE 


SPEED 


VARIABLE  DESCmPTION 


wind  speed  at  the  surface 


UNITS/ 


meters/sec 


LNRHL 


SHEAR  I  wind  shear  between  high  and  low  layers 


low  level  average  temperature 


middle  level  average  temperature 


temperature  difference  between  sigma  layers  1  and  0 
(approx.  Tgjo  ^  —  Tsuriace) 


temperature  difference  between  sigma  layers  4  and  1 
(approx.  T500  mb  ~  Tgso  mb) 


low  level  average  relative  humidity 


natural  logarithm  of  low  level  relative  humk^ty 


tow  level  relative  humidity  squared 


tow  level  relative  humidity  raised  to  the  power  of  4 


middle  level  average  relative  humidity 


LNRHMO  I  natural  logarithm  of  middle  level  relative  humidity 


middle  level  relative  humidity  squared 


middle  level  relative  humidity  raised  to  the  power  of  4 


maximum  relative  humidity  of  six  sigma  layers 


relative  hurr‘dity  at  layer  above  maximum 


3-hr  change  in  low  level  relative  humidity 


3-hr  change  in  middle  level  relative  humidity 


G-hr  change  in  low  level  relative  humidity 


G-hr  change  in  middle  level  relative  humkfity 


9-hr  change  in  low  level  relative  humidity 


G-hr  change  in  middle  level  relative  humidity 


mb  thickness 


surface  pressure 


3-hr  change  in  surface  pressure 


G-hr  chap  ^  •*  surface  pressure 


9-hr  change  in  surface  pressure 


model  surface  height 


low  level  divergence 


meters/sec 


D3RHL 


03RHM 


D6RHL 


D6RHM 


D9RHL 


09RHM 


[0.0.1.01 


[0.0,1.01 


[0.0,1.01 


[0.0,1.01 


[-1.0,1.01 


[-1.0,1.01 


[-1.0.1.01 


[-1.0,1.01 


[-1.0,1.01 


[-1.0,1.01 


Pascals 


Pascals 


Pa'A^' 


Pascals 


second"’ 
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VARIABLE 

NAME 

VARIABLE  DESCRIPTION 

UNITS/ RANGE 

OVM 

middle  level  divergence 

second"’ 

OVH 

high  level  divergerK» 

second"’ 

D3DVL 

3-hr  change  in  low  level  divergence 

second"’ 

D3DVM 

3-hr  change  in  middle  level  divergefx:e 

secorxt"’ 

D3DVH 

3-hr  change  in  high  level  divergence 

second"’ 

D6DVL 

6-hr  change  in  low  level  divergence 

second"’ 

D6DVM 

6-hr  change  in  middle  level  divergence 

second"’ 

06DVH 

G-hr  change  in  high  level  divergence 

second"’ 

D9DVL 

9-hr  change  in  low  level  divergence 

second"’ 

DSOVM 

9-hr  change  in  middle  level  divergence 

SeCOTKt"’ 

D9DVH 

9-hr  change  in  high  level  divergence 

secoTKl"’ 

VTL 

low  level  vorticity 

secortd"’ 

VTM 

middle  level  vorticity 

SeCOTKj"’ 

VTH 

high  level  vorticity 

second"’ 

D3VTL 

3-hr  change  in  low  level  vorticity 

second"’ 

D3VTM 

3-hr  change  in  middle  level  vorticity 

second"’ 

D3VTH 

3-hr  change  in  high  level  vorticity 

second"’ 

D6VTL 

6-hr  change  in  low  level  vorticity 

second"’ 

D6VTM 

6-hr  change  in  middle  level  vorticity 

secorxJ"’ 

D6VTH 

6-hr  change  in  high  level  vorticity 

second"’ 

D9VTL 

9-hr  change  in  low  level  vorticity 

second"’ 

D9VTM 

9-hr  change  in  middle  level  vorticity 

secorxt"’ 

D9VTH 

9-hr  change  in  high  level  vorticity 

second"’ 

TAL 

low  level  temperature  advection 

Kelvin/sec 

TAM 

middle  level  temperature  advection 

Kelvin/sec 

TAH 

high  level  temperature  advection 

Kelvin/sec 

D3TAL 

3-hr  change  in  low  level  temperature  advection 

Kelvin/sec 

D3TAM 

3-hr  change  in  middle  level  temperature  advection 

Kelvin/sec 

D3TAH 

3-hr  change  in  high  level  temperature  advection 

Kelvin/sec 

D6TAL 

6-hr  change  in  low  level  temperature  advection 

Kelvin/sec 

D6TAM 

6-hr  change  in  middle  level  temperature  advection 

Kelvin/sec 
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VARIABUE 

NAME 

VARUBLE  DESCRIPTION 

D6TAH 

6-hr  change  in  high  level  temperature  advection 

Kelvirr/sec 

D9TAL 

9-hr  change  in  low  level  temperature  advection 

Kelvin/sec 

D9TAM 

9-hr  change  in  middle  level  temperature  advection 

Keivin/sec 

D9TAH 

9-hr  change  in  high  level  temperature  advection 

Keivin/sec 

WL 

low  level  vertical  velocity 

secorxl"’ 

WM 

middle  level  vertical  velocity 

second"’ 

WH 

high  level  vertical  velocity 

second*’ 

D3WL 

3-hr  change  in  low  level  vertical  velocity 

second"’ 

D3WM 

3-hr  change  in  middle  level  vertical  velocity 

second"’ 

D3WH 

3-hr  change  in  high  level  vertical  velocity 

second"’ 

D6WL 

6-hr  change  in  low  level  vertical  velocity 

secorrd"’ 

D6WM 

6-hr  change  in  middle  level  vertical  velocity 

second"’ 

D6WH 

6-hr  change  in  hi^  level  vertical  velocity 

secorxJ"’ 

09WL 

9-hr  change  in  low  level  vertical  velocity 

second"’ 

09WM 

9-hr  change  in  middle  level  vertical  velocity 

secorrd"’ 

D9WH 

9hr  change  in  high  level  vertical  velocity 

second"’ 

WRHL 

product  of  low  level  vertical  velocity  and  relative  humidity 

second"’ 

WRHM 

product  of  middle  level  vertical  velocity  and  relative  humidity 

second"’ 

DELEV 

elevation  gradient  in  local  average  wind  direction  (normalized  to 
grid  box  size) 

dimensionless 

ACCL3 

3-hr  upwind  NEPH  cloud  cover  (averaged  to  half  mesh)  using  average 
low  level  wind 

[0.0,1. 0) 

ACCL6 

6-hr  upwind  NEPH  cloud  cover  (averaged  to  half  mesh)  using  average 
low  level  wind 

[0.0,1.01 

ACCL9 

9-hr  upwind  NEPH  cloud  cover  (averaged  to  half  mesh)  using  average 
low  level  wind 

[0.0,1.01 

ACCM3 

3-hr  upwind  NEPH  cloud  cover  (averaged  to  half  mesh)  using  average 
middle  level  wind 

[0.0,1.01 

ACCM6 

6-hr  upwind  NEPH  doud  cover  (averaged  to  half  mesh)  using  average 
middle  level  wind 

[0.0,1.01 

ACCM9 

9-hr  upwind  NEPH  cloud  cover  (averaged  to  half  mesh)  using  average 
middle  level  wind 

[0.0,1.01 
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VARIABLE 

NAME 

VARIABLE  DESCRIPTION 

UNITS/  RANGE 

ACCH3 

3-hr  upwind  NEPH  doud  cover  (averaged  to  half  mesh)  using  average  high 
level  wind 

(0.0,1.01 

ACCH6 

6-hr  upwind  NEPH  cloud  cover  (averaged  to  half  mesh)  using  average  hi^ 
level  wind 

(0.0,1.01 

ACCH9 

9-hr  upwind  NEPH  cloud  cover  (averaged  to  halt  mesh)  using  average  high 
level  wind 

(00,1.0] 

RH50 

binary  variable  that  is  1  if  relative  humidity  is  greater  than  50%,  else  0 

0,1 

RH70 

binary  variable  that  is  1  if  relative  humidity  is  greater  than  70%,  else  0 

0,1 

RHQO 

binary  variable  that  is  1  if  relative  humidity  is  greater  than  90%,  else  0 

0,1 

NEPH3 

RTNEPH  cloud  cover  at  the  eighth-mesh  grid  point  at  the  center  (lower 
right)  of  the  half-mesh  cell  3  hours  before  valid  time 

(0.0,1.01 

NEPH6 

RTNEPH  cloud  cover  at  the  eighth-mesh  grid  point  at  the  center  (lower 
right)  of  the  half-mesh  cell  6  hours  before  valid  time 

(0.0,1.01 

NEPH9 

RTNEPH  cloud  cover  at  the  eighth-mesh  grid  point  at  the  center  (lower 
right)  of  the  half-mesh  ceil  9  hours  before  valid  time 

(0.0,1.01 

NEPH12 

RTNEPH  cloud  cover  at  the  eighth-mesh  grid  point  at  the  center  (lower 
right)  of  the  half-mesh  cell  12  hours  before  valid  time 

(0.0,1.01 

NEPH15 

RTNEPH  cloud  cover  at  the  eightti-mesh  grid  point  at  the  center  (lower 
right)  of  the  half-mesh  cell  15  hours  before  valid  time 

(0.0,1.01 

NEPH18 

RTNEPH  cloud  cover  at  the  eighth-mesh  grid  point  at  the  center  (lower 
right)  of  the  half-mesh  ceil  18  hours  before  valid  time 

(0.0,1.01 

NEPH21 

RTNEPH  cloud  cover  at  the  eighth-mesh  grid  point  at  the  center  (lower 
right)  of  the  half-mesh  cell  21  hours  before  valid  time 

[0.0,1.01 

NEPH24 

RTNEPH  cloud  cover  at  the  eighth-mesi;  grid  point  at  the  center  (lower 
right)  of  the  half-mesh  cell  24  hours  before  valid  time 

(0.0,1.01 
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APPENDIX  B 


In  this  appendix  we  present  a  derivation  of  the  formula  for  posterior  probability 
that  produced  the  evaluation  results  on  MOS  forecast  probabilities  reported  in  Sec¬ 
tion  3.5.3.  While  elementary,  the  derivation  requires  careful  use  of  delta  functions  to  rep¬ 
resent  probability  measures  (equivalently,  probability  density  functions)  which  are 
defined  on  a  continuous  domain  (the  interval  [0, 1]  in  this  case)  yet  have  non-zero  weight 
(or  probability)  at  an  isolated  point. 

In  the  derivation,  we  assume  that  the  event  of  interest,  E(i,  j),  has  been  fixed  (i.e., 
a  ccf  bin,  Bi,  and  a  threshold  value,  Tj,  have  been  selected)  and  suppress  notation  identify¬ 
ing  these  parameters.  Define 

p  =  the  true  probability  that  E(i,j)  will  occur 
p  =  the  MOS  forecast  probability  of  E(i,j)  based  on  data  set  Dj 

n  =  the  number  of  observations  in  data  set  Dg  for  which  p  is  in  bin  Bj 
k  =  the  number  of  observations  in  data  set  Dg  for  which  E(i,  j)  occurs . 


Now,  our  Bayesian  viewpoint  picks,  somewhat  arbitrarily,  the  following  prior  prob¬ 
ability  measure  on  the  unknown  (true)  probability,  p: 


flp) 


0.5  6(p-p)  +  0.5,  0  ^  p  s  1.0 

0  ,  otherwise 


(B-1) 


where  5(x)  is  a  delta  function  concentrated  at  x  =  0.  Once  data  set  D2  has  been  analyzed 
(i.e.,  the  values  of  n  and  k  have  been  determined  for  D2),  the  likelihood  of  a  particular 
value,  y,  is  defined  by 

Uy/D^)  =  y‘‘(l-y)"-^  (B-2) 


Note  that  the  likelihood  has  the  form  of  a  binomial  density  function  (lacking  only  a  scaling 
constant)  because  data  set  D2  is  assumed  to  provide  n  independent  trials.  Each  “trial”  con¬ 
sists  of  comparing  the  actual  cloud  cover  fraction  at  the  forecast  valid  time  with  the 
selected  threshold,  Tj. 
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Using  Eqs.  B-1  and  B-2,  the  posterior  probability  measure  on  p  has  the  form 


^p/d/P^  - 


(0.5c)  p‘Hl-p)"-‘^[6(p-p)  +  1.0],  0  5  p  ^  1.0 
0  ,  otherwise 


(B-3) 


where  the  constant,  c,  must  be  chosen  so  that  the  integral  over  [0, 1]  is  equal  to  1.0.  Carry¬ 
ing  out  the  integration  results  in 


1 


0.5c  =1  p*‘(l-p)"  ‘‘  -h  I  x^Hl-xl^'^^dx. 


J  x*Hl-x)"‘*‘dx.J 


(B.4) 


So,  to  calculate  the  posterior  probability  (p*)  that  p  =  p,  we  can  evaluate  the  following 

p+A 

p*  =  lim  I  fp/D/P)dp  •  (B-5) 

A  — 0 


p-A 


Substituting  Eqs.  B-3  and  B-4  into  Eq.  B-5  and  evaluating  the  limit  yields  the  desired  for¬ 
mula,  Eq.  3-25. 
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